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Nomenclature 


Speed  of  sound 

Local  coefficient  of  friction 
Specific  heat  at  constant  pressure 
Velocity  ratio, 

e 2 

Enthalpies,  defined  In  the  expression  H = h + 

Thermal  conductivity 

Eddy  conductivity 

Characteristic  problem  dimension,  length  of  the  model  in 
question 

Defined  In  Eq  (19-1) 

Mach  number 
Prandtl  number 
Pressure 

Heat  flux  or  heat  flow  per  unit  area 
Gas  constant,  1716  ft2/sec2R  for  air 
Reynolds  number 

Radial  coordinate  (body  radius)  for  the  case  of  the  axl- 
symmetrlc  cone,  measured  perpendicularly  from  the  longitudinal 
centerline,  fig  3 

Viscosity  constant  of  Sutherland  (198.6  R) 

Nondimenslonal  position,  x/L 

Stanton  number,  ■■■  — — jr-r 
peue(He_V 

Temperature 

Transverse  curvature  term  equal  to  ■— 

Veloclty  component  along  (perpendicular  to)  the  streamwlse 
direction 

Transformed  velocity  expression  defined  In  Eq  (18-3) 

Defined  In  Eq  (27) 


x,y  Body  surface  oriented  coordinate  system  in  which  x runs 

parallel  to  the  stream  direction,  pointing  downstream,  and 
y is  perpendicular  to  x and  is  directed  into  the  external 
flow 

Greek  Symbols 

a Defined  in  Eq  (19-2) 

3 Defined  in  Eq  (19-3) 

r Streamwise  Intermittency  distribution  or  probability  factor 

Y The  gas  constant,  ratio  of  specific  heats 

Y*  The  intermittency  factor  of  Klebanoff 

A Change  in  variable  quantity 

6 Boundary  layer  thickness 

6*  Displacement  thickness 

e Eddy  viscosity 

e"  Eddy  viscosity  function  defined  following  Eq  (22) 

A 

e Eddy  viscosity  function  defined  following  Eq  (22) 

(pv)w 

X A nondlmens lonal  mass  transfer  rate,  t— -t 

^pu'»  or  e 

x\  Transformed  perpendicular  boundary  layer  coordinate  and  non- 

dimensional  distance  along  this  coordinate 

9 Static  temperature  ratio,  J- 

e 

9 Momentum  thickness 

p Molecular  viscosity 

v Kinematic  viscosity, 

P 

5 Transformed  streamwise  boundary  layer  coordinate  and  nondl- 

mensional  length  along  this  coordinate 

p Density 

t Shear  Stress 

u Exponent  of  the  viscosity  law  of  Sutherland 
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Subscripts  and  Superscripts 

V 

e Condition  at  the  edge  of  the  boundary  layer,  also  indicative 

of  the  input  or  environmental  conditions  for  Itract  in  the 
cone  study 

® Free  stream  or  unperturbed  condition 

j Flow  index,  j = 1 for  conical  flow,  j = 0 for  flow  over  a 

flat  surface 

6*(9)  When  used  with  Re,  denotes  Reynolds  number  based  on  displace- 
ment thickness  (momentum  thickness) 

o Total  or  stagnation  condition  except  for  rQ 

Primed  quantities  indicate  instantaneous  departures  from  a 
mean  state  or  condition  in  the  turbulence  model.  The  accom- 
panying bars  over  the  primed  symbols  denote  a time  averaged 
quantity. 

ref  Reference 

t Turbulent  condition 

w Condition  at  the  surface  of  the  plate  or  cone 

x Denotes  a particular  real  x station  along  the  surface  of  the 

model 
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Abstract 

This  study  followed  the  work  of  Or  J.  Shang,  Flight  Dynamics  Labor- 
atory, Wright  Patterson  Air  Force  Base,  Ohio.  Given  a Fortran  code 
written  by  Dr  Shang  that  solved  for  the  characteristics  of  a Laminar, 
transitional,  and  turbulent  boundary  layer,  the  problem  was  to  modify 
the  existing  program  to  predict  the  boundary  layer  over  a flat  plate  and 
sharp  nosed  axisymmetric  cone  with  mass  transfer  as  a boundary  condition 
at  the  surface  of  the  model.  The  surface  of  the  model  was  maintained  at 
a constant  temperature,  and  only  the  cases  in  which  air  was  the  trans- 
ferred gas  were  studied. 

To  solve  this  problem  the  boundary  layer  was  described  by  the 
standard  boundary  layer  equations  for  continuity,  momentum,  and  energy. 
Incorporating  mass  transfer  as  a boundary  condition,  the  governing  equa- 
tions underwent  the  transformation  of  Probstein-Elliott  and  Levy-Lees. 

The  resulting  equations  and  boundary  conditions  were  solved  by  finite 
differencing  techniques  for  nondimensionalized  velocity  components  and 
temperature  at  a finite  number  of  nodes  in  the  boundary  layer  field  of  flow. 

To  verify  the  modified  code,  three  studies  were  performed.  First, 
the  code  was  verified  using  analytical  and  some  experimental  data  from 
Schllchtlng  for  laminar,  subsonic  flow  over  a flat  surface  with  constant 
suction.  Second,  the  code  was  verified  for  turbulent,  subsonic  flow  over 
a flat  surface  with  constant  suction  to  the  asymptotic  suction  limit  and 
for  small  rates  of  blowing,  using  experimental  results  from  Moffat  and 
Kays.  Finally,  the  code  was  verified  with  mixed  success  for  hypersonic 
laminar,  transitional,  and  turbulent  flow  over  an  axisymmetric  cone  for 
small  rates  of  blowing  using  the  experimental  results  of  Martelluccl, 
Laganelli.and  Hahn. 
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The  surface  of  the  model  was  maintained  at  a constant  temperature,  and  only 
the  cases  In  which  air  was  the  transferred  gas  were  studied. 
t To  solve  this  problem  the  boundary  layer  was  described  by  the  standard 
boundary  layer  equations  for  continuity,  momentum,  and  energy.  Incorporating 
mass  transfer  as  a boundary  condition,  the  governing  equations  underwent  the 
transformation  of  Probsteln-Elllott  and  Levy-Lees.  The  resulting  equations 
and  boundary  conditions  were  solved  by  finite  differencing  techniques  for 
nondlmenslonallzed  velocity  components  and  temperature  at  a finite  number  of 
nodes  In  the  boundary  layer  field  of  flow.<^- 

To  verify  the  modified  code,  three  studies  were  performed.  First,  the 
code  was  verified  using  analytical  and  some  experimental  data  from  Schllchtlng 
for  laminar,  subsonic  flow  over  a flat  surface  with  constant  suction.  Second, 
the  code  was  verified  for  turbulent,  subsonic  flow  over  a flat  surface  with 
constant  suction  to  the  asymptotic  suction  limit  and  for  small  rates  of 
blowing,  using  experimental  results  from  Moffat  and  Kays.  Finally,  the  code 
was  verified  with  mixed  success  for  hypersonic  laminar,  transitional,  and 
turbulent  flow  over  an  axl symmetric  cone  for  small  rates  of  blowing  using 
the  experimental  results  of  Martelluccl,  Laganelll,  and  Hahn. 
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Calculating  the  effects  on  boundary  layer  flows  subjected  to  mass 
transfer  perpendicular  to  the  surface  has  provided  engineers  a rela- 
tively Inexpensive  model  to  study  ablative  effects.  This  model  has 
provided  a means  by  which  to  study  the  heating  effects  at  the  surface, 
skin  friction,  and  the  effects  on  the  boundary  layer  profiles.  The 
purpose  of  this  paper  was  to  Investigate  the  effects  of  this  mass 
transfer  at  the  solid  boundary  by  means  of  a numerical  code. 

Definition  of  the  Problem  and  Purpose  of  the  Study 

The  Flight  Dynamics  Laboratory  (FDL)  possessed  a digital  computer 
code,  called  Itract,  which  computed  the  characteristics  of  laminar  and 
turbulent  boundary  layers  over  flat  plates  and  axi symmetric,  conical 
bodies  for  the  case  with  no  mass  transfer  at  the  surface.  To  Initiate 
this  computation  the  following  quantities  were  specified  as  Inputs: 
gamma,  the  gas  constant;  the  Prandtl  numbers,  both  laminar  and  turbulent; 
free  stream  mach  number,  static  temperature,  and  density;  the  exponent 
of  the  viscosity  law  of  Sutherland;  a temperature  ratio,  wall  tempera- 
ture to  free  stream  stagnation  temperature;  a point  of  transition  from 
laminar  to  turbulent  flow  along  the  surface;  and  a flagged  quantity 
which  specified  eddy  model  zero  or  eddy  model  one  for  computation  of 
the  eddy  viscosity.  With  these  Inputs,  Itract  provided  a description 
of  boundary  layer  features.  Some  of  the  output  of  Interest  In  this 
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study  Included  the  following:  the  local  mach  number  for  any  point  in 

the  field,  boundary  layer  thickness,  displacement  thickness,  momentum 
thickness,  the  coefficient  of  friction,  eddy  viscosity,  a Stanton  number 
descriptive  of  heat  transfer  at  the  surface,  and  boundary  layer  profiles 
for  velocity,  static  temperature,  and  density. 

Starting  with  the  original  code  of  FDL  this  study  was  divided  into 
three  sequential  steps.  The  first  step  was  to  learn  as  much  about  the 
computer  code  as  possible.  This  step  Included  a study  of  the  key 
equations  of  motion,  energy,  and  continuity  needed  for  boundary  layer 
solution.  The  second  step  was  to  incorporate  the  needed  changes  into 
the  code  that  would  include  the  new  boundary  condition  of  mass  transfer 
at  the  surface  of  the  body  exposed  to  an  environment  of  fluid  flow.  The 
last  step  was  to  verify  the  change  by  comparing  key  output  predictions 
of  the  computer  code  with  the  results  of  analytical  expressions  presented 
in  Schlichting  and  laboratory  experiments  for  studies  of  flow  over  a 
flat  plate  and  flow  over  a slender,  axlsymmetrlc  cone.  Completing  these 
three  steps,  the  purpose  of  this  study  was  to  extend  the  usefulness  of 
a turbulent  boundary  layer  code  by  Incorporating  a change  that  would 
allow  consideration  of  mass  transfer  as  a boundary  condition,  and 
thereby,  study  Its  effect  on  boundary  layer  characteristics. 

Scope  and  Assumptions 

In  defining  the  area  of  study  the  topic  was  limited  and  the 
following  assumptions  were  made.  First,  boundary  layer  computations 
and  comparisons  were  performed  on  flat  surfaces  and  axlsymmetrlc  cones 
with  sharp  leading  edges  or  tips.  For  both  models  there  were  negligible 
effects  due  to  the  stagnation  region  at  the  leading  edges,  and  In  the 
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case  of  the  plate,  the  shocking  phenomena  was  neglected.  Shapiro 
alluded  to  the  validity  of  this  assumption  of  free  stream  conditions 
existing  some  distance  downstream  of  the  leading  edge  of  a plate  In 
fig  28-21(c)  and  subsequent  text  (Ref  1:1149-1150).  Eckert  discussed 
this  Idea  further  as  mach  numbers  reached  supersonic  and  higher  (Ref  2: 
10-11).  Thus,  free  stream  conditions  were  assumed  to  exist  downstream 
of  the  shock  wave.  Further,  the  angle  of  incidence  of  the  models  was 
assumed  to  be  zero  with  respect  to  the  flow  in  the  free  stream.  Second, 
consideration  was  given  only  to  the  cases  of  air,  at  surface  temperature, 
being  blown  through  the  surface  into  the  boundary  layer,  or  the  boundary 
layer  of  air  flow  being  sucked  through  the  porous  surface  into  the 
model.  Thus,  this  study  did  not  Include  the  effects  of  chemical  reactions 
that  might  occur  by  mixing  nonsimilar  gases  In  the  boundary  layer. 

Taking  the  transferred  gas  to  be  at  the  temperature  of  the  wall,  which 
was  assumed  constant,  helped  to  limit  and  simplify  the  problem  and  the 
transfer  model.  This  was  a realistic  limit  as  many  experiments  in  wind 
tunnels  were  performed  under  these  constant  temperature  conditions. 

Third,  only  small  rates  of  Injection  or  suction  were  compared  with 
experimental  results,  although  the  limiting  transfer  rates  of  the  code 
were  Investigated.  It  was  assumed  that  mass  transfer  effects  were  con- 
fined to  the  boundary  layer  (Ref  3:1, 5-6).  The  solution  of  this 
problem  was  based  partially  on  the  boundary  layer  equation  of  motion  of 
Prandtl.  To  have  considered  massive  transfer  rates  would  have  violated 
the  proposition  of  Prandtl  that  6 was  much  less  than  the  characteristic 
length  of  the  model.  Thus,  a fourth  stipulation  was  that  6 would  be 
much  less  than  L.  Further,  the  pressure  change  across  this  boundary 
layer  thickness  was  negllblble,  and  considered  zero  In  the  analytical 
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solution.  Fifth,  the  problem  was  limited  to  experimental  cases  where 
r pressure  change  along  the  stream  was  also  negligible.  This  was  con- 

sistent with  the  two  models  studied.  Numerically,  dp/dx  was  considered 
zero.  Finally,  in  the  studies  of  both  the  flat  plate  and  the  conical 
flow,  the  mass  transfer  rate  was  considered  constant  over  the  region 
of  transfer  unless  indicated  otherwise  in  the  experiment.  Also,  the 
flow  was  considered  fully  turbulent  throughout  the  length  of  the  model 
unless  another  transition  point  was  clearly  indicated  in  the  experi- 
mental results.  This  list  of  items  provided  the  limits  and  scope  of 
this  study.  The  following  chapter  presents  a background  of  information 
relevant  to  this  study. 

i o 


II.  A Background  of  Information 


Interest  In  boundary  layers  perturbed  by  mass  transfer  at  the  sur- 
face has  been  evident  from  numerous  laboratory  experiments  In  which  a 
model  equipped  with  a surface  blowing  apparatus  was  exposed  to  the  free 
stream  environment  of  a wind  tunnel.  More  recently,  computer  codes  have 
been  designed  to  compute  the  same  fluid  characteristics  as  measured  in 
the  experimental  efforts.  In  both  these  studies  those  features  of  the 
boundary  layer  that  were  of  greatest  interest  included  the  following: 

a)  Boundary  layer  velocity  profile  shape, 

b)  Energy  (temperature)  profile  shape, 

c)  Thicknesses  - boundary  layer,  displacement,  and  momentum, 

d)  Skin  friction  reduction  for  the  blowing  case,  and 

e)  Heat  transfer  blockage  for  the  blowing  case. 

In  the  experimental  study  these  features  have  been  obtained  by 
measuring  a restricted  number  of  quantities. 

The  devices  used  to  measure  these  quantities  in  experiments  on 
boundary  layers  have  Included  heat  transfer  gages,  pressure  sensors, 
temperature  probes,  and  mass  Injection  concentration  probes  (Ref  4: 
1-10,  32-35,  46-51).  The  same  quantities  measured  by  these  devices 
have  been  computed  by  analytical  methods.  Such  a method  or  computer 
code  was  written  for  the  Flight  Dynamics  Laboratory,  Wright  Patterson 
Air  Force  Base,  Ohio. 

This  code  was  written  to  obtain  numerical  solutions  of  the 
governing  turbulent  boundary  layer  equations.  Because  of  limited 
understanding  of  turbulent  processes,  completely  general  solutions  to 
these  equations  have  not  been  possible.  By  use  of  an  empirical  eddy 
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viscosity  model  of  these  processes,  however,  the  system  of  governing 
equations  was  solved  directly.  The  basic  eddy  model  used  in  this  study 
was  that  of  Cebeci,  Smith,  and  Mosinskis.  The  model  assumed  an  inner 
and  an  outer  viscous  layer  within  the  boundary  layer.  The  expression 
for  e In  the  inner  region  was  based  on  the  mixing  length  theory  of 
Prandtl  as  follows: 


•i  - * It 


O) 


where  l was  equal  to  K^y.  To  account  for  the  region  close  to  the  wall. 
Van  Driest  offered  a modification  to  the  mixing  length  of  Prandtl.  The 
new  expression  for  t was 


L « Kjy  (1  - exp(-y/A))  (2) 

where  A was  equal  to  26v  (T,/Pw)”^2*  The  shear  stress  close  to  the  wall 
was  written 

(t)y  (3) 

If  A were  redefined  to  26v  (t/p)"1^2,  then  expanded 


A * 26v 


t . J-l/2 
JL  + 2E  Z, 

P dx  p 


Finally  then. 


2..2 


’inner  ” V 


1 - exp 


‘sfcr 


— + 

P dx  p 


1/2' 


$ 


The  expression  for  e In  the  outer  region  was 

h C (Vu><<y 


“outer 


(4) 


(6) 
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This  became  eddy  model  zero  In  the  code  and  differed  from  eddy  model 
one  which  was  formed  by  multiplying  eouter  or  eddy  model  zero  by  the 
intermittency  factor  of  Klebanoff, 


For  6 defined  as  the  distance  from  the  surface  to  a point  in  the  field 
at  which  u was  equal  to  .995U*,,  studies  have  shown  that  the  value  for 
K<|  was  .4  and  the  value  for  Kg  was  .0168  (Ref  5:1975-1976).  Having 
selected  the  model  of  Cebecl  to  describe  turbulent  activity  within  th£ 
boundary  layer,  there  remained  the  problem  of  solving  the  system  of 
governing  equations  for  laminar,  transitional,  and  turbulent  compressible 
boundary  layers  (Ref  6). 

Finite  differencing  techniques  were  incorporated  to  obtain  solu- 
tions of  the  governing  system  for  both  flat  plate  and  axi symmetric 
conical  flows.  The  numerical  technique  Involved  a simultaneous  solution 
of  the  equations  of  momentum,  energy,  and  continuity  by  a tridiagonal 
matrix  inversion  routine.  Through  an  iterative  procedure,  the  solutions 
of  all  three  were  brought  into  convergent  harmony  yielding  results 


which,  otherwise,  would  have  been  gained  only  through  laboratory  experi- 
ments. Some  of  the  mathematical  modeling  Incorporated  with  these  three 
governing  equations  Included  a two-layer  concept  within  the  turbulent 
boundary  layer  with  appropriate  eddy  viscosity  models  used  for  the  Inner 
and  outer  regions.  These  models  were  considered  in  addition  to  the 
molecular  viscosity  term  applicable  In  laminar  flow.  Further,  a speci- 
fied turbulent  Prandtl  number  related  turbulent  heat  flux  to  the  Reynolds 
stress.  Finally,  mean  properties  within  the  transition  region  between 
laminar  and  turbulent  flow  were  computed  by  multiplying  the  eddy 
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viscosity  by  an  Intermlttency  factor  that  characterized  the  growth  rate 
or  production  of  turbulence  within  a flow  whose  origin  was  laminar  (Ref 
7;  Ref  8:1-4).  With  these  models  Incorporated,  the  solution  followed. 

The  Equations  to  b£  Solved 

The  flow  of  compressible,  viscous,  heat  conducting  fluid  was 
described  by  the  equations  of  continuity,  Navier-Stokes,  and  energy, 
together  with  a supporting  equation  of  state,  a heat  conductivity  law, 
and  the  viscosity  law  of  Sutherland.  To  arrive  at  such  a description 
was  to  accept  the  propositions  of  Prandtl.  Osborne  Reynolds  was  the 
first  to  study  turbulent  flow  in  1883.  He  said  that  the  Instantaneous 
fluid  velocity  satisfied  the  Navier-Stokes  equations,  and  that  this 
velocity  was  comprised  of  a mean  velocity  and  a fluctuating  component. 
He  modified  the  Navier-Stokes  expressions  with  these  fluctuating 
components,  called  Reynolds  Stresses,  and  by  making  boundary  layer 
approximations  he  presented  the  governing  equations  as  follows  (Ref  8: 
11-12): 

Continuity 

4 (rJpu)  * 57 

Momentum 
Energy 

p[“&(cpT)  * (v  ♦ v <cpt)]  ■ u 3?  * 7 k & (CPT) 

HSsP’Hp  (-cpp7Tr>]-pr7ri7  (io> 


d£  . 1 3 J f 

•a7+^1y  r l 


u!y  + 


pu'v' 


(9) 
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Appendix  A was  Included  for  further  clarification  of  the  above  system 
(Ref  9:145-150). 

Eqs  (8),  (9),  and  (10),  valid  descriptions  for  laminar  and  turbu- 
lent flow,  were  the  laminar  governing  equations  with  the  addition  of 
turbulent  fluctuating  quantities  which  represented  the  apparent  turbu- 
lent mass,  shear,  and  heat  flux  terms.  These  turbulent  additions  were 
Incorporated,  again,  through  mathematical  modeling.  The  apparent  mass 
flux,  p'v',  was  incorporated  by  the  new  variable,  v;  the  apparent  shear 
stress,  pu ' v ' , became  part  of  the  eddy  model;  and  the  apparent  heat  flux, 
CpPV 1 t 1 , was  modeled  through  an  eddy  conductivity  term,  Ky.  These 
relationships  were  defined  by  the  following  equations: 


v = v + 


e = - p 


3u/3y 


v * T ' 

kT  = ~cp  p FT73y 

With  these  quantities  incorporated,  the  perfect  gas  law  and  the  viscosity 
relation  of  Sutherland  were  also  added: 


Perfect  gas  law 


P - cr 


Viscosity  law 


Te  + S 

ttt 


(air  only) 


Thus,  the  system  of  governing  equations  to  be  solved  consisted  of  three 
nonlinear  partial  differential  equations  and  two  algebraic  expressions. 
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( 


But  In  the  present  form  this  system  had  a singularity  at  x equal  to  zero, 
the  leading  edge.  To  alleviate  this  singularity,  and  to  reduce  the 
growth  of  the  boundary  layer  as  the  solution  proceeded  downstream  for 
numerical  efficiency,  a variable  transformation  was  made  (Ref  8:13-15). 

The  Transformed  Plane 

The  transformation  of  Probstein-Elliott  and  Levy-Lees  was  used  in 
this  analytic  study.  The  transformation  was  written  as  follows: 

ZM  - f PeVero2jdx  (14) 


p u r ^ fV  . 

n(x,y)  - fe  jr°-  t j -£-  dy  (15) 

Jo  pe 

Next,  the  relation  between  derivatives  in  the  real  (x,y)  plane  and  the 
transformed  plane  (5,n)  followed: 


p u u r^ 
e e e o 


peuerojtj  mm 


(16) 

(17) 


Then,  the  three  parameters,  F,  0,  and  V were  defined  as  follows: 


F 


£ 

V 


JL 


p u y r ^ 
We  o 


PvrpjtJ 

m.  t 


(18) 
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With  this,  the  final  working  form  of  the  governing  system,  prior  to 
linearization,  was  reached.  Further  definitions  included 


l = _£H_ 


pepe 


'P  e 


6 - 

P ' ue  dS 


(19) 


Finally,  the  solvable  form  of  the  governing  system  was  obtained  as 
follows: 

Continuity 

li+2«t?+F*0  (20) 

Momentum 

2?Flitvl^-H(t23®l)  + 6<F2-i>'0  (21) 


Energy 


Casting  Eqs  (20),  (21),  and  (22)  into  a finite  difference  form,  this 
system  represented  a means  by  which  a boundary  layer  could  be  studied 
numerically.  With  the  Inclusion  of  boundary  conditions,  this  system 
was  solvable.  For  purposes  of  this  study  the  boundary  conditions  were 
as  follows:  F(?,0)  » 0 

VU.O)  « vw(0 

9(£,0)  ■ 0W,  a constant 
F(C. ne)  - ! 

0(E,ne)  - 1 (Ref  8:13-18;  Ref  10) 


(23) 


This  chapter  has  Introduced  the  boundary  layer  problem,  and  methods 
by  which  this  problem  has  been  rtudied  and  solved.  The  methods  presented, 
experimental  and  analytical,  represented  the  techniques  employed  by  those 
in  the  engineering  community  who  have  studied  boundary  layer  flow  exten- 
sively. The  numerical  solution  ultimately  depended  on  the  boundary 
conditions  Imposed  on  the  differential  equations.  Further,  the  boundary 
condition,  Eq  (23-2),  was  to  become  the  primary  area  of  study  for  this 
thesis.  This  quantity,  VWU),  would  ultimately  provide  Itract  with  the 
capability  to  investigate  the  effects  of  mass  transfer  on  a boundary 
layer.  The  original  FDL  code  solved  the  boundary  layer  problem  for  no 
mass  transfer,  or  VwU)  equal  to  zero.  With  a Vw(£)  model  Incorporated 
to  simulate  the  mass  transfer  of  air,  the  code  would  solve  the  boundary 
layer  problem  such  as  that  investigated  by  the  experimental  study  men- 
tioned at  the  beginning  of  this  chapter.  To  better  understand  this 
numerical  solution  It  was  necessary  to  include  a program  description. 
Chapter  III. 


The  computer  code,  I tract,  solved  the  system  of  nonlinear  para- 


bolic partial  differential  equations,  Eqs  (20),  (21),  and  (22),  by 
casting  this  system  into  a series  of  linear  finite  difference  expres- 
sions. Coincidentally,  the  transformation  from  the  real  (x,y)  plane  to 
the  (£.n)  plane  cast  the  boundary  layer  Into  a rectangular  grid  of  nodes 
with  the  surface  of  the  model  located  at  the  level  j = l,  as  shown  in 
fig  1. 


Fig.  1.  Finite  Difference  Grid  for  Boundary  Layer  (From  Ref  8:33) 

The  solution  of  this  system  of  finite  difference  equations  was  approxi- 
mated by  computing  values  of  F,  9,  and  V at  each  of  the  nodes  within 
the  grid.  With  values  for  these  variables  at  stations  1-2  and  1-1  the 
values  of  F,  0,  and  V were  solved  at  station  1 from  the  surface  to  the 
edge  of  the  boundary  layer  using  a three-point  differencing  scheme  and 
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a tridiagonal  matrix  Inversion  routine.  With  the  boundary  layer  solu- 
tion completed  at  station  i the  problem  was  stepped  in  the  streamwise 
direction,  £,  to  station  i+1  and  the  node  by  node  computation  was  per- 
formed again  from  the  surface  to  the  edge  of  the  boundary  layer.  The 
entire  program  was,  therefore,  a sequential  solution  of  a series  of 
columns  of  nodes  from  the  leading  edge  to  the  trailing  edge  of  the  sur- 
face or  model.  For  the  particular  problem  considered  in  this  study, 
the  program  followed  the  step-by-step  procedure  depicted  in  fig  2,  with 
a program  listing  included  in  Appendix  B. 


Inputs  including  Y,  Pr,'  Pr*,  M^,  Tm,  stepping  data 
In  transformed  plane,  ui,  Tv,/T0,  turbulent  transition 
and  intermittency  consideration,  L,  p^,  and  mass 
transfer  conditions 

1 ■ 

| Compute  Reoo  | 


I 

Continued  on  following 
page 


15 


Fig.  2.  Flow  Diagram  of  the  Logical  Steps  to  Solution  within  Itract 


r 


4 


I 


Key  portions  of  the  foregoing  logic  required  further  explanation.  There- 
fore, Appendix  C was  Included  to  discuss  four  Important  subsystems  of 
the  original  code.  These  subsystems  included  nondimensionalization  of 
the  working  variables  and  Initialization  of  the  grid,  generation  of  the 
finite  difference  system,  the  computation  of  e,  and  the  compution  of 

Ste  and  Cf  . A Fortran  computer  code  key  was  also  Included  in  Appendix 
e 

D.  With  an  understanding  of  these  features  of  the  code,  the  boundary 
condition  of  mass  transfer  was  considered.  Including  this  boundary  con- 
dition represented  the  major  modification  to  the  original  code,  and  the 
remainder  of  this  chapter  was  devoted  to  an  explanation  of  this  addition. 


Mass  Transfer 


Mass  transfer  at  the  surface  was  defined  by  the  expression  (pv)w. 


Consistent  with  the  nondimensionallzed  variables  used  in  this  problem  a 

(pv)w 

mass  transfer  factor,  j --y - , was  defined  and  used  to  express  the 

'p  • or  e 

amounts  of  mass  transfer  being  considered  in  any  particular  problem.  This 
transfer  model  was  incorporated  through  the  variable  transformation 


JL 
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Wero  - 


pvr  jt^| 
F|^+— 2- 
3x 


(18-3) 


and  expressed  In  the  equation  of  continuity 


(20) 


where,  V appeared  explicitly  In  the  finite  difference  expression  for 
continuity. 

Considering  Eq  (18-3)  in  detail,  the  following  points  were  noted: 
First,  at  the  surface  F or  u/uft  was  zero.  Second,  t^,  where  t was  the 
ratio  r/r0,  was  set  equal  to  one.  This  assumption  was  made  following 
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the  proposition  that  <5  was  much  less  than  the  radius  of  the  cone.  Figure 
3,  though  the  boundary  layer  was  shown  out  of  proportion,  depicted  the 
pictorial  justification  for  this  assumption. 


% D 1 V 1 

Third,  from  an  earlier  definition  restated,  v was  equal  to  v + . 

It  was  noted  that  vw  was  equal  to  vw  at  the  wall  or  surface  as  the  appar- 
ent mass  flux,  p ' v ' , was  zero.  With  these  three  propositions  Eq  (18-3) 
was  expressed  for  the  wall  condition  as  follows: 

(pv) 

V r = V(5,0)  or  V(i,l)  in  the  grid  notation  (24) 

w p u u r J 
Me  eMe  o 

Returning  to  the  entire  equation  of  continuity,  integration  yielded  an 
expression  for  V at  any  grid  point  at  station  s^: 

V(1.j)  ■ V(1,l)  - (2S  ||+F)1  dn  (25) 

where  V (1 ,1 ) was  the  boundary  condition  of  mass  transfer.  To  Include 
or  V(1,l),  further  substitution  within  Eq  (24)  was  performed.  From 
fig  3,  rQ  was  equal  to  x sin  $c>  and  It  could  be  shown  that 
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and  with  one  additional  intermediate  step  it  was  shown  that 


and,  finally, 
Vw 


p u V L J-SI 

^oo  oo  oo  y 


1/2 


Lj(sin<j>c)j  (2X)1/2 


im 


1/2 


pe  ue  ye 
Ip»  u«o  Pref 


P v 
w w 


P<oUJ 


I 

(26) 


(27) 


(28) 


(29) 


Now  in  terms  of  quantities  immediately  available  in  the  program,  this 
expression  was  cast  into  an  equivalent  form  using  nondimensional  program 
variables  (Ref  8:18,35;  Ref  10).  With  Eq  (29)  including  the  effects  of 
mass  transfer,  the  equation  of  continuity  was  considered  next. 

Cast  Into  a form  of  finite  differences,  continuity  was  expressed  as 
follows: 


C3(l ,3)V(1 ,j+l ) + C3(l , 2 ) V ( 1 ,j)  + C3(l ,1 )V(1 ,j-l ) + A3(l ,2)F(1,j ) = D3(l) 

(30) 


At  the  surface  this  expression  simplified  to 


C3(1.1)V(1,J-1)  « D3(l ) 


(31) 


Setting  C3(1,1)  equal  to  one  and  D3(1)  equal  to  the  right  side  of  Eq 
(29)  the  mass  transfer  boundary  condition  had  been  set  and  was  Included 
with  the  other  boundary  conditions  In  solving  the  system  of  finite 
difference  equations. 
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In  order  to  set  an  appropriate  boundary  condition  at  each  station 
along  5 during  the  computation,  two  subroutines  were  added  to  the  pro- 
gram. For  the  case  in  which  a constant  mass  transfer  rate  was  specified 
in  a real  sense  along  the  surface  from  some  Initial  longitudinal  station 
to  a second  station  where  mass  transfer  was  terminated,  subroutine 
Conblw  provided  an  appropriate  transformed  value  for  the  transfer  at 
each  station  computed.  A second  subroutine,  Genblw,  provided  the  same 
information,  but  for  a generally  varying  mass  transfer  rate.  Using  a 
linear  Interpolation  between  stations  of  known  mass  transferring  strength, 
the  boundary  condition  was  computed  for  each  streamwise  station  within 
the  specified  region  of  mass  transfer.  Finally,  although  not  incorpor- 
ated Into  the  code,  an  approximation  using  a cubic  spline  description 
between  known  or  specified  points  of  transfer  rate  was  devised  during 
this  study.  It  was  thought  that  this  technique  would  have  provided  a 
better  description  of  a generally  varying  mass  transfer  rate, and  the 
theory  of  the  proposed  modification  was  included  in  Appendix  E (Ref  11). 
However,  with  the  other  modifications  completed,  numerical  solutions 
with  mass  transfer  were  compared  with  analytical  and  actual  experimental 
results,  and  the  results  of  those  comparisons  were  included  in  Chapter 
IV. 
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IV.  Results  and  Discussion  of  the  Study,  Flat  Plate  and  Cone 

1 

The  modified  program  was  compared  with  theory  and  data  from  three 
primary  sources.  First,  using  mostly  analytical  expressions  and  some 
experimental  data  presented  in  Schlichting  (Ref  12)  a study  was  made  of 
laminar,  subsonic  flow  over  a flat  plate  for  the  cases  of  no  mass  trans- 
fer and  a constant  rate  of  suction  throughout  the  length  of  the  model. 
Second,  from  the  results  of  an  experiment  performed  by  Moffat  and  Kays 
(Ref  13)  a comparison  was  made  for  fully  turbulent,  subsonic  flow  over 
a flat  plate.  The  comparison  was  made  for  the  cases  of  no  mass  transfer, 
constant  blowing,  and  constant  suction  over  a specified  region  of  the 
model.  Finally,  from  an  experiment  performed  for  flow  over  a sharp  nosed, 
axisymmetric  cone  by  Martellucci,  Laganelll,  and  Hahn  (Ref  14;  Ref  15), 
data  was  obtained  to  test  the  computer  code  for  the  case  of  hypersonic 
flow.  For  this  case  of  hypersonic,  conical  flow,  the  numerical  results 
were  compared  in  laminar,  transitional,  and  turbulent  environments  for 
the  cases  of  no  mass  transfer  and  positive  mass  transfer  or  blowing. 

In  these  studies  a number  of  important  assumptions  were  made,  some 

of  which  were  mentioned  earlier  In  Introductory  comments.  The  boundary 

layer  thickness,  6,  was  minutely  small  compared  to  the  characteristic 

length,  L.  The  velocity  gradient,  |y,  was  large  in  this  region,  and 
du 

the  shear  stress,  assumed  large  values.  Beyond  the  boundary  layer 
no  large  velocity  gradients  existed  and  viscosity  was  negligible.  The 
flow  was  considered  Invlscld  and  potential  beyond  the  edge  of  the 
boundary  layer.  Finally,  the  Navler-Stokes  equations  were  simplified 
to  the  boundary  layer  equations  to  describe  flow  characteristics  for  y 
( less  than  «(Ref  12:117-121). 
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Schlichtlng,  Primarily  an  Analytical  Verification  for  Laminar  Flow  Over 
a Flat  Plate 

For  purposes  of  this  study,  a hypothetical  model  and  some  flow  con- 
ditions were  needed  to  make  the  comparison  between  analytic  results  and 
the  predictions  of  the  code.  A comparison  for  the  case  of  no  mass  trans- 
fer was  followed  by  a study  with  a constant  rate  of  suction  over  a flat 
plate. 

Beginning  with  the  case  of  laminar  subsonic  flow  with  no  mass 
transfer  at  the  surface,  working  variables  were  assigned  the  following 
values.  Re  was  adjusted  to  about  1.(10)^  in  keeping  with  the  laminar 
propositions  of  Blaslus.  Further,  Tw  was  selected  equal  to  to  be 
consistent  with  the  environment  for  which  Eq  (33)  would  be  valid.  It 
was  also  consistent  with  the  results  of  Eq  (32); 

u2ro  _ 

Tw  * T«  < 'R*  : Heat  wal1  «-  9”  (32) 

The  right  side  of  this  inequality  for  the  test  under  Investigation  pro- 
duced an  extremely  small  difference  between  and  T^,  and  hence,  there 
was  zero  heat  transfer  or  the  adiabatic  case.  Finally,  a length  of  three 
ft  was  chosen  for  the  hypothetical  model  of  the  flat  plate  In  order  to 
specify  Rew.  The  remaining  Inputs  for  this  first  test  for  program  veri- 
fication Included  an-M^  equal  to  .01,  a T^  of  533.1  R,  and  a pw  equal 
» lbfsec* 

to  1.12(10)”*  — —~i — . For  verification  In  at  least  this  case  of  steady 
ft4 

laminar  flow  over  a flat  plate  without  mass  transfer,  the  resulting 
computations  at  station  s equal  to  .155  and  station  s equal  to  .750  were 
chosen  for  comparison  with  the  calculations  of  the  exact  expressions 
, listed  In  Schlichtlng.  The  quantities  chosen  for  comparison  were  6*, 

u 

t , cf,  and  6*/0.  Also  Included  was  a comparison  of  velocity 
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and  temperature  profiles  with  data  presented  in  Schlichting  from  the 
work  of  Hantzsche  and  Wendt  (Ref  12:323).  F -om  Schlichting,  the 


following  expressions  of  Blasius  were  used  for  computation: 


1.721  = 6* 


.332 


Ux) 


vx 


'vx' 


VU0 


'f  - w'"'  - .332  I v 


pu; 


1/2 

1/2 
1/2 


u«x. 


.332 
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(Rex) 
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The  results  of  a comparison  between  computations  performed  by  the 
use  of  the  above  expressions  and  by  calculations  performed  by  the  com- 
puter code,  Itract,  were  summarized  in  Table  I. 


Table  I 

A Comparison  of  Methods  for  Boundary  Layer  Calculations 


Quantity  for 
Comparison 

Station  s = .155 

Station  s = .750 

■ 

Schlichting 

Itract 

Schlichting 

Itract 

6*(ft) 

g 

|| 

SI 

SI 

wimv 

EQQBil 

5.52  (10)"4 

5.52  (10)'4 

cf 

B3DSSI 

SI 

—a 

7.71  (10)"4 

6*/9 

2.59 

2.61 

2.59 

2.62 

To  quantify  the  difference  noted  between  the  predictions  of  Itract  and 
the  analytical  or  experimental  data,  an  error  was  defined  as  the  quotient 
of  the  absolute  difference  between  the  quantities  compared  and  the  larger 


of  the  two  quantities.  Thus,  the  results  of  Table  I demonstrated  a close- 
ness to  within  the  following  percentage  errors.  At  station  s equal  to 
.155  the  calculations  of  6*  were  within  1.5  percent,  t results  were 
nearly  identical,  the  calculations  of  c^  were  within  1.1  percent,  and 
the  computations  of  6*/G  were  within  .8  percent  of  one  another.  A similar 
trend  was  noted  at  station  s equal  to  .750.  The  calculations  of  6*  were 

within  1.5  percent,  t..  results  were  again  equivalent,  the  calculations 

w 

of  Cf  were  within  .6  percent,  and  the  computations  for  6*/Q  were  within 
1.1  percent  of  one  another. 

Further  tests  for  verification  of  the  program  in  this  first  case 
study  were  accomplished  by  comparing  velocity  and  temperature  profiles 
calculated  by  Hantzsche  and  Wendt  with  the  predictions  of  Itract  (Ref 
12:323,  fig  13-11).  It  was  noted  that  5(y/<5)  in  the  code  was  equivalent 
to  the  n of  Blasius.  Further,  the  — of  Blasius  was  equivalent  to 

Uoo 

.995  -p  in  Itract.  With  these  relationships  plus  the  computational 
equivalence  of  Tg  in  Itract  to  T^  in  Schlichting,  the  results  of  the  com- 
parison were  listed  in  Table  II  with  a graphical  presentation  of  the 
velocity  profiles  presented  in  fig  4.  Concerning  the  velocity  profile, 
the  data  of  station  s equal  to  .731  was  used  for  comparison,  but  with 
similarity  of  solution  for  this  particular  Investigation  and  the  non- 
dimensional  Ized  nature  of  the  data,  another  station  would  have  been 
equally  valid  for  comparison. 
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Table  II 


A Comparison  of  Velocity  and  Temperature  Profiles 


Blasius 

n 

Itract 

y/6 

H and  W 

U/U  co 

Itract 
.995  u/ue 

H and  W 
T/T 

00 

Itract 

T/Te 

1 

.2 

.35 

.36 

1.0 

1.0 

2 

.4 

.64 

.65 

1.0 

1.0 

3 

.6 

.84 

.85 

1.0 

1.0 

4 

.8 

.95 

.96 

1.0 

1.0 

5 

1.0 

.99+ 

.99+ 

1.0 

1.0 

The  greatest  error  in  this  comparison  was  less  than  2.8  percent  within 
the  velocity  profile  study.  With  these  profile  comparisons  the  investi- 
gation for  the  first  case  was  completed.  Case  two  added  mass  transfer 
to  the  problem. 

Initial  testing  of  the  actual  modification  to  the  program  began  with 

the  addition  of  a small  mass  transfer  condition,  constant  suction.  Kays 

also  presented  the  method  of  Rubesin  for  analytically  studying  large 

mass  transfer  rates  (Ref  16:324-325).  To  complete  the  study  for  small 

constant  suction  the  experimental  and  analytical  work  of  Head  and  Iglisch, 

as  published  in  Schllchting,  was  used  to  verify  the  results  of  Itract 

(Ref  12:373,  Fig  14.11.1).  T remained  equal  to  T for  this  second  test. 

pw  was  then  shown  to  be  equal  to  pm  by  the  equation  of  state,  and  from 

(pvlw  , ,-4 

fig  14.11.1,  therefore,  was  equal  to  -1.6(10)  in  Itract.  Data 

was  collected  at  the  nondimensional  streamwise  position 


( 


.077 

“fid 


-v. 


wl 


(34) 
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This  implied  that  the  profile  data  of  Head  was  recorded  along  the  flat 
surface  at  a station  where  Rex  was  approximately  3 . 00 ( 1 0 )6 . For  this 
comparison,  then,  the  hypothetical  length  of  the  model  was  extended  from 
3 ft  to  30  ft,  where  Reynolds  numbers  of  this  size  would  be  encountered. 
Laminar  conditions  were  still  assumed  to  exist.  Assuming  in  fig  14.11.1 
of  Schlichting  that  6 was  approximately  1.8  mm,  a graphical  comparison 
for  this  test  was  presented  in  fig  5. 

To  show  the  effect  on  the  shape  of  the  velocity  profile  by  the  addi- 
tion of  suction,  fig  6 portrayed  the  results  of  Itract  for  the  boundary 
layer  flows  with  and  without  suction.  These  results  agreed  with  the 
results  presented  in  Schlichting  (Ref  12:369,  fig  14.6). 

Now,  as  with  the  first  case  study,  there  existed  an  exact  solution 
for  flow  over  a flat  plate  with  continuous,  constant  suction.  The 
following  equation  represented  an  exact  solution  of  the  complete  Navier- 


Stokes  equations: 


V 

u(y)  ■ 1 - exp  — - 


for  v(x,y)  = vw<  0 


From  this  expression  came  the  following  equations: 


Tw  " P(-Vw>u«’  and  hence. 


T -2v 
w _ w 


f l/2pu  l uo 
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It  was  noted  that  in  each  case,  g*,  0,  and  c^  were  constant.  This  solu- 
tion was  realized  only  at  some  distance  from  the  leading  edge.  The 
boundary  layer  grew  from  zero  at  the  leading  edge  and  continued  down- 
stream asymptotically  to  the  values  predicted  by  Eq  (36).  These  values 
were  reached  at  what  was  termed  the  asymptotic  suction  layer  limit. 
Iglisch  has  shown  that  the  asymptotic  state  was  reached  after  a length 
of  about 


x = 


4v 


-v 


(37) 


w 


To  simulate  this  asymptotic  solution  the  length  of  the  hypothetical 
model  was  extended  still  further  to  3000  ft,  and  the  remaining  input  con- 
ditions were  held  constant.  Iglisch  then  predicted  an  asymptotic  solution 
by  station  s equal  to  .156.  Itract  had  come  within  2.3  percent  of  the 
final  asymptote  by  s equal  to  .155.  Table  III  summarized  the  results  from 
the  equations  of  the  exact  solution  above,  and  compared  those  calculations 
with  the  corresponding  predictions  of  Itract  at  an  s of  .347,  the  point  of 
closest  approach  to  the  analytical  asymptotic  values. 


Table  III 

A Laminar  Flat  Plate  Study  with  Suction 


( 


Quantity 

Exact 

Solution 

Itract 

Percent 

Error 

6*/L 

6.25(10)'* 

6.04(10)'® 

3.3 

0/L 

3.12(10)'6 

2.97(10)'6 

4.9 

cf 

3.20(10)"4 

3.24(10)'4 

1.2 
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Finally,  all  testing  thus  far  that  included  mass  transfer  had  been 
accomplished  using  the  routine  that  incorporated  constant  mass  transfer 
rates  at  the  surface.  Before  investigating  other  experiments  with  flat 
plates,  the  variable  mass  transfer  routine  was  verified.  First,  using 
the  three  ft  model,  Itract  computed  a boundary  layer  perturbed  by  a con- 
stant rate  of  suction  from  a point  one  ft  from  the  leading  edge  to  a point 
two  ft  from  the  leading  edge.  The  computation  was  repeated  with  the  same 
Inputs  with  the  exception  that  the  variable  mass  transfer  routine  was 
called  to  compute  the  boundary  condition  in  lieu  of  the  constant  mass 
transfer  routine.  Identical  results  were  noted  for  the  two  tests. 

With  this  final  check  the  verification  process  departed  from  the 
laminar  flow  study  and  considered  turbulent  flow  over  a flat  surface. 

For  this  study  the  results  of  experiments  performed  by  Moffat  and  Kays 
were  used. 

Moffat  and  Kays , A Verification  for  Turbulent  Flow  Over  a.  Flat  Plate 
Using  Experimental  Results 

R.  J.  Moffat  and  W.  M.  Kays  of  Stanford  University  performed  an 

experiment  In  which  they  were  primarily  concerned  with  heat  transfer 

through  a turbulent  boundary  layer  over  a flat  plate  which  was  perturbed 

by  both  positive  and  negative  mass  transfer  at  the  solid  boundary.  The 

results  of  their  wind  tunnel  study  provided  a criterion  for  evaluating 

the  heat  transfer  model  of  Itract  under  turbulent  conditions.  Heat 

transfer  In  the  experiment  was  quantified  In  the  form  of  a Stanton  number, 

St.  The  accuracy  of  the  apparatus  used  allowed  determination  of  the 

-4 

Stanton  numbers  to  within  10  units  over  most  of  the  range  of  mass 
transfer.  The  experiment  was  performed  on  a transfer  range  from  the 
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asymptotic  suction  layer  limit  discussed  earlier  to  the  apparent  blow  off 
or  separation  of  the  boundary  layer.  Presented  in  this  section  are 
results  of  testing  and  a discussion  of  a parameter  study  performed  to 
minimize  the  effects  of  higher  order  terms  not  included  and,  hopefully, 
match  this  numerical  model  with  the  experimental  environment  for  the  no 
transfer  case.  With  accurate  predictions  for  this  case,  the  results  for 
small  amounts  of  blowing  and  suction  were  given  next.  Finally,  the  range 
of  accurate  prediction  of  the  computer  code  was  tested,  with  these  results 
included  last. 

To  begin,  a wind  tunnel  run  was  chosen  with  the  following  conditions: 
u^  was  equal  to  44.5  ft/sec,  was  524.0  R,  and  Tw  was  556.6  R.  The 

experimental  data  collected  was  listed  in  Table  V.  The  length  of  the 
model  was  given  as  8 ft.  It  was  assumed  that  the  last  value  of  Rex  was 
taken  from  the  end  of  the  plate,  and  could  bs  considered  a close  approxi- 
mation to  Re  . Further,  it  was  assumed  that  the  flow  was  turbulent  over 

CD 

the  entire  length  of  the  wind  tunnel  model.  A parameter  study  was  then 
begun  to  find  the  best  combination  of  those  variables  which  described  the 
grid  to  minimize  error  caused  by  the  truncation  of  higher  order  terms,  and 
pick  two  parameters  which  helped  describe  the  characteristics  of  the  flow. 
These  two  classes  of  variables  included  the  following:  XXK,  the  constant 

ratio  of  any  two  successive  An  spaces;  PRT,  the  turbulent  Prandtl  number 
taken  to  be  1.  or  .9  in  the  literature;  XINTER  set  to  1.  or  0.  depending 
on  whether  eddy  model  one  or  eddy  model  zero  was  to  be  used;  DYW,  the  size 
chosen  for  An^ ; and  IED6E,  the  total  number  of  divisions  in  n to  be  used  in 
the  computation  of  the  grid.  The  objective  was  to  closely  predict  the 
Stanton  number  for  a corresponding  Rex  that  ranged  from  4.55(10)*,  where 
measurements  of  heat  transfer  began,  to  the  end  of  the  plate  at  an  Rex  of 
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2.14(10)®.  Table  IV  of  Appendix  F summarized  the  combinations  of  vari- 
ables with  Table  V of  that  same  appendix  actually  presenting  the  results 
of  those  variable  combinations.  The  figures  of  column  4 produced  the 

best  match  with  the  experimental  results.  Excluding  the  readings  at  an 

4 5 

Rex  of  4.55(10)  the  greatest  error  was  recorded  at  an  Rex  of  2.27(10) 

with  an  error  of  5.7  percent.  Column  3 had  produced  nearly  identical 

results,  but  had  incorporated  Inefficiently  small  stepping  increments 

Into  the  numerical  scheme.  A graphical  presentation  of  the  experimental 

results  with  the  analytical  predictions  of  column  4 and  column  1 was 

Included  In  fig  7.  In  a final  note,  with  the  exception  of  readings  at 

Rex  values  of  4.55(10)4,  2.27(10)5,  and  3.18(10)5,  the  remaining  errors 

were  less  than  or  equal  to  3.9  percent. 

With  the  case  for  zero  mass  transfer  recorded,  two  more  experimental 

runs  were  Investigated.  First,  an  experiment  which  Included  a blowing 

(pv)  _3 

rate,  |~uy,  1 - (10)  was  run  under  the  following  conditions:  u^  was 

equal  to  44.1  ft/sec,  T^  was  525.7  R,  and  Tw  was  557.7  R.  In  a simula- 
tion by  Itract  the  results  were  presented  In  Table  VI  of  Appendix  F with 
a graphical  presentation  Included  In  fig  8.  From  Table  VI  It  was  noted 
that  In  setting  XINTER  equal  to  0.,  and  thereby  using  eddy  model  zero  In 

the  calculation  of  e,  more  accurate  Stanton  numbers  resulted.  Next,  an 

(pv)w  ,3 

experiment  which  Included  a rate  of  suction,  ^puy  , equal  to  -1.15(10) 
was  run  under  the  following  conditions:  u^  was  equal  to  42.5  ft/sec,  T^ 
was  524.3  R,  and  T^  was  549.7  R.  Again,  the  results  of  a simulation  by 
Itract  were  presented  In  Table  VII  of  Appendix  F with  the  graphical 
equivalent  Included  In  fig  8.  Unlike  the  case  with  blowing  the  tabular 
results  for  this  case  with  suction  showed  that  the  more  accurate  pre- 
dictions of  Stanton  numbers  came  by  setting  XINTER  equal  to  1.,  and 


Fig.  7.  Matching  Itract  with  Moffat  and  Kays  for  the  Case  of 
Zero  Mass  Transfer  (Ref  13) 
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Fig.  8.  A Comparison  of  I tract  and  Moffat  and 
Kays  for  the  Cases  of  Blowing  and  Suction  (Ref  13) 


thereby  using  eddy  model  one.  Finally,  with  simulations  performed  for 

both  the  small  positive  and  negative  mass  transfer  cases,  it  was  then 

appropriate  to  find  the  limits  of  accurate  simulation  by  Itract. 

In  this  final  phase  of  flat  plate  testing  Itract  was  simulated  at 

the  extreme  limits  of  the  Moffat  and  Kays  experiment.  In  the  limiting 

case  for  suction,  termed  the  asymptotic  suction  layer  limit,  Itract  was 

able  to  predict  Stanton  numbers  to  within  5.3  percent,  excluding  one 

reading  taken  at  a station  where  Rex  was  equal  to  4.3(10)*.  The  wind 

tunnel  conditions  for  this  test  Included  the  following:  u was  equal  to 

00 

41.8  ft/sec,  T was  523.8  R,  and  T..  was  552.  R.  The  rate  of  suction  was 


w -3 

equal  to  -7.77(10)  . The  tabular  results  of  this  test  were 


(pv) 

(pu)_ 

included  in  Table  VIII  with  the  graphical  summary  included  in  fig  8. 

Again,  as  with  lower  suction  rates,  more  accurate  results  were  noted 

when  eddy  model  one  was  used.  However,  unlike  the  case  for  suction,  in 

the  testing  of  positive  mass  transfer  or  blowing,  Itract  was  unable  to 

predict  heat  transfer  to  the  limiting  point  of  blow  off  or  boundary  layer 

(pv)  _3 

separation,  which  occurred  experimentally  near  ^uy-  equal  to  9.6(10)  . 

The  results  of  the  predictions  of  Itract  for  rates  of  blowing  equal  to 
l.(lO)"'*  have  already  been  presented.  For  the  code,  the  limiting  trans- 
fer  rate  for  which  there  existed  experimental  data  was  1.91(10)  . At 

this  transfer  rate  the  numerical  scheme  could  compute  the  boundary  layer 
problem  without  an  error  finish.  The  results  of  this  test  were  included 
In  Table  IX  with  a graphical  summary  Included  In  fig  8.  It  was  noted 
that  with  the  finer  mesh  of  nodal  points  Itract  was  able  to  predict 
consistently  the  Stanton  numbers  for  various  Rex  up  to  a point  where  the 
numerical  scheme  failed.  While  the  scheme  was  able  to  compute,  Itract 
consistently  predicted  Stanton  numbers  3.(10)"*  less  than  the  experimental 
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results  from  an  Rex  of  2.28(10)5  to  1.23(10)6  where  the  program  exper- 
ienced an  error  finish.  With  a coarser  mesh  of  nodal  points  Itract  was 
able  to  complete  the  numerical  computation,  but  with  predictions  of 
Stanton  number  that  were  not  as  close  as  previous  tests.  Rather  than  a 
nearly  constant  difference  of  prediction  as  previously  seen,  the  results 
of  this  test  showed  Itract  to  predict  Stanton  numbers  lower  than  experi- 
mental by  about  22.8  percent.  From  an  Rex  of  2.28(10)^  through  the  end 
of  the  computation  the  greatest  deviation  from  this  figure  was  to  25.7 

percent.  Finally,  in  a test  case  for  a mass  transfer  of  3.8(107  , using 

Ari-j+i  _4 

a coarse  grid  of  ^ equal  to  1.15,  An-j  at  5.(10)  , and  100  divisions 

3 

In  the  n grid,  Itract  was  able  to  successfully  compute  the  boundary 
layer  without  error  finish.  However,  experimental  values  of  Stanton 
number  ranged  from  2.36(10)"^  to  6.2(10)"^,  and  with  Itract  predicting 
values  consistingly  greater  than  5.(10)"^  below  the  experiment,  the 
results  were  not  included. 

The  results  for  blowing  equal  to  1.(10)  displayed  the  limit  of 
positive  mass  transfer  rate  with  which  Itract  could  compute  accurately. 

.3 

Beyond  a transfer  rate  of  3.8(10)  Itract  was  neither  able  to  predict 
Stanton  numbers  nor  successfully  complete  the  computations  without  an 
error  finish.  This  completed  the  comparison  with  the  experiment  by 
Moffat  and  Kays. 

Martel! uccl , Laganelll , and  Hahn,  A Study  of  Turbulent  Flow  Over  an 
Axl symmetric  Cone  with  Experimental  Results 

A.  Martel 1 uccl,  A.  L.  Laganelll,  and  J.  Hahn  of  the  General  Elec- 
tric Reentry  and  Environmental  Systems  Division  performed  an  experiment 
over  a two  year  period  In  which  they  were  concerned  with  heat  transfer 
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behavior  and  boundary  layer  profile  characteristics  for  hypersonic  flow 
over  a sharp  nosed,  slender,  axi symmetric  cone.  Their  experimental 
results  of  heat  transfer  and  profile  data  provided  numerous  quantities 
by  which  to  evaluate  the  modified  code. 

In  the  experiment,  data  was  collected  for  nominal,  positive  mass 
transfer  rates  as  follows:  0.,  5.(10)  1.(10)  and  1.5(10)"^.  All 

four  transfer  rates  were  investigated  in  this  study,  with  comparisons 
between  data  and  numerical  predictions  made  for  the  heat  transfer  at 
the  surface,  the  velocity  profile,  and  the  static  temperature  profile. 

In  making  this  comparison  there  was  a problem  in  describing  the  flow 
environment  downstream  of  the  leading  oblique  shock  wave. 

Unlike  the  study  of  flow  over  a flat  plate,  the  oblique  shocking 
effect  on  the  cone  was  great  enough  to  significantly  change  the  fluid 
state  downstream  of  the  shock  wave.  Therefore,  for  purposes  of  computa- 
tion, the  actual  free  stream  conditions  were  not  of  direct  use  to  the 
computer  code.  Rather,  the  environment  downstream  of  the  shock  wave  was 
the  needed  condition  for  input  Into  Itract.  Computing  these  conditions 
for  Input  would  have  been  a time  consuming  problem  in  itself,  and  the 
needed  additions  to  the  existing  code  to  perform  this  computation  were 
not  pursued.  In  order  to  provide  the  conditions  at  the  edge  of  the 
boundary  layer,  graphs  of  characteristics  of  flow  over  a cone,  such  as 
those  found  in  NACA  1135,  were  considered.  Not  only  did  the  resolution 
of  the  graphical  Information  seem  Inadequate  for  the  range  of  mach 
number  being  considered,  but  the  data  presented  was  for  an  Invlscid, 
compressible  solution.  Tabulated  data  such  as  that  Included  In  reference 
17  was  considered,  and  though  accurate.  It  still  posed  data  for  an 
Invlscid  solution.  Investigations  were  made  using  the  data  of  the 


Inviscld  solution  in  reference  17  as  inputs  to  Itract.  It  was  judged 
that  this  method  did  not  yield  results  close  enough  to  the  physical 
situation  at  hand  to  be  considered  a valid  approximation.  To  obtain 
viscous  inputs  for  Itract,  the  decision  was  made  to  use  data  presented 
with  the  results  of  the  experiment  performed  at  General  Electric. 

A review  of  the  experimental  technique  was  appropriate.  As  stated 
previously,  the  data  collected  in  the  experiment  was  of  two  categories. 
These  two  categories  of  data  were  collected  in  separate  runs  of  the  wind 
tunnel.  Initially,  the  model  of  the  cone  was  exposed  to  flow  at  an 
equal  to  eight  for  a few  seconds.  The  heat  transfer  data  was  collected 
and  flow  within  the  tunnel  was  stopped.  After  the  surface  data  had  been 
taken,  flow,  again  at  an  equal  to  eight,  was  started.  The  inter- 

action of  the  flow  over  the  model  of  the  cone  was  allowed  to  reach  an 
equilibrium  state,  and  the  second  category  of  data  profile  Information 
was  collected  (Ref  4:11).  Within  this  profile  data,  the  following  mea- 
surements or  computations  were  taken  for  various  stations  along  the  cone 
Me,  Tg,  ue,  Pg,  (pv)w,  (Pu)e,  and  Tw/TQ.  The  above  quantities,  mostly 
representative  of  conditions  at  the  edge  of  the  boundary  layer,  became 
the  new  conditions  at  infinity  to  be  used  as  inputs  to  Itract.  These 
Inputs  were  used  by  Itract  to  predict  surface  as  well  as  the  field  data 
of  the  boundary  layer.  With  this  assumption,  the  following  approxima- 
tions were  made  for  computational  purposes:  First,  where  data  from 

multiple  stations,  both  longitudinal  and  azimuthal,  along  the  model  was 
catalogued  for  the  same  wind  tunnel  environment,  an  arithmetic  average 
of  quantities  such  as  Me,  Tg,  and  pg  at  these  stations  was  used  to  com- 
pute a new,  constant  M^,  T^,  and  p^  for  Itract.  Further  It  was  approxi- 
mated that  T was  a constant  ratio  equal  to  an  arithmbtic  average  of 
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the  readings  taken  along  the  surface  in  a streamwise  direction.  In  fact, 
wall  temperature  did  vary  in  the  experiment  and  the  temperature  ratio 
was  seen  to  vary  plus  or  minus  three  or  four  percent  from  the  figure  used 
in  computation.  It  was  noted  that  one  term  in  the  denominator  of  the 
expression  used  to  compute  Stanton  numbers  was  (1  -T  /T  ),  and  values  for 
Tw/T0  of  .5  to  .8  were  conmon  (Ref  14;  Ref  15).  Also,  since  the  defin- 
ition of  the  Stanton  number  of  Martellucci  was  actually  an  St  , it  was 

co 

P U 

necessary  to  multiply  the  Itract  figure  by  the  factor  -pp  prior  to  com- 
parison  with  the  experimental  data.  Finally,  there  were  three  descriptions 
for  mass  transfer  rate:  First,  a nominal  figure  for  blowing  was  presented 

such  as  5.(10)"^,  1.(1Q)”*,  and  1.5(10)”'*.  Second,  an  actual  measurement 
of  this  blowing  rate  would  be  found  by  performing  the  division  (pv).,/(pu)  . 

W 00 

This  was  designated  as  A . In  like  manner,  (pv)  /(pu)  was  computed  and 

W C 

defined  as  Ag.  All  three  had  different  actual  values,  and  all  three 

figures  were  tested  in  the  modified  code.  Though  all  were  describing  the 

same  mass  transfer  activity,  A0  was  finally  selected  as  the  appropriate 

boundary  condition  for  this  code. 

Using  the  assumptions  and  approximations  listed  above,  the  cases 

tested  and  presented  were  of  four  categories:  First,  a study  of  the 

case  for  no  mass  transfer  was  considered.  After  this,  three  investiga- 

tions  followed  with  nominal  mass  transfer  rates  of  1.5(10)  , 1.(10)  , 

»-4 

and  5.(10)  . These  four  cases  comprised  the  entire  study  of  flow  over 

the  sharp  nosed,  axisynmetrlc  cone. 

Beginning  the  study  of  flow  over  a cone  with  a nonblowing  case,  an 
experimental  test  case,  data  group  132,  was  chosen  from  the  results  of 
Martellucci,  Laganelll,  and  Hahn.  This  was  a data  group  depicting  heat 
transfer  at  the  surface  of  the  cone  In  the  form  of  St  for  numerous 

oo 


longitudinal  positions  along  the  surface.  Connected  with  this  heat  trans- 
fer data  group  were  data  groups  74  through  79  that  presented  profile  or 
field  data  and  were  the  product  of  the  same  free  stream  conditions  as 

data  group  132.  The  free  stream  conditions  included  an  M equal  to  7.87, 

_5  lbfsec^ 

T equal  to  about  92.9  R,  and  a p of  2.59  (10)  j — . Using  data 

fC 

groups  74  through  79,  the  actual  Input  conditions  to  Itract  were  an  Mg 

equal  to  6.84,  a T equal  to  121.  R,  and  a equal  to  3.66(1 0)-5 
lbfsec2  e e 

t — . The  length  of  the  model  was  five  ft,  the  point  of  transition 

ft* 

was  approximated  from  the  experimental  Stanton  number  curve  to  be  about 
1.33  ft  from  the  tip  of  the  cone,  and  an  average  Tw/TQ  was  found  to  be 
approximately  .68.  Using  this  Information  a tabulated  comparison  of 
the  heat  transfer  results  was  listed  in  Table  X of  Appendix  G with  a 
graphical  depiction  Included  in  fig  9.  This  graph  not  only  showed  the 
results  of  Itract  in  comparison  with  the  experimental  data  but  provided 
theoretical  boundaries  for  heat  transfer  as  predicted  by  Bell  Aircraft 
Corporation  (Ref  18).  The  lower  Bell  curve  predicted  heat  transfer 
assuming  the  flow  was  laminar  throughout  the  length  of  the  model.  The 
upper  Bell  curve  predicted  the  heat  transfer  assuming  fully  turbulent 
flow  for  the  entire  length  of  the  model.  Concerning  the  prediction  of 
Itract,  It  was  noted  that  the  curve  continually  overpredicted  the  experi- 
mental heat  transfer,  followed  similar  heat  transfer  trends  as  the  flow 
proceeded  along  the  surface,  and  settled  to  within  2.3  to  8.4  percent 
of  the  data  for  the  last  1.5  ft  of  the  cone.  It  was  further  found  that, 
unlike  the  flat  plate  study  with  blowing,  eddy  model  one  yielded  the 
better  results  In  predicting  heat  transfer  for  the  cone.  Some  of  the 
disparity  of  heat  transfer  prediction  In  the  region  of  transition  was 
due  to  an  approximated  turbulent  transition  point.  The  first  departure 
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Station  (in.) 

g.  9.  Matching  Itract  with  Martel! ucci,  Laganelli,  and  Hahn  for 
the  Case  of  Zero  Mass  Transfer,  Heat  Transfer  Group  132  (Ref  14) 


from  a linear  trend  in  the  Stanton  number  data  plotted  on  a logarithmic 
scale  was  used  as  the  point  of  transition  (Ref  10).  To  further  inves- 
tigate this  case  for  no  mass  transfer  two  more  cases  were  considered. 

It  was  thought  at  General  Electric  that  the  results  of  two  partic- 
ular cases  offered  excellent  references  or  test  cases  by  which  to 
compare  the  predictions  of  Itract  (Ref  19).  The  first  case  was  tenned 
a hot  wall  experiment,  a nearly  adiabatic  wall,  and  was  similar  to  each 
of  the  succeeding  cases  with  mass  transfer  that  would  be  studied.  The 

free  stream  conditions  for  this  test,  data  group  150,  included  an  Mro  of 

_c  lbfsec2 

8.0,  a T^  of  97.6  R,  and  a of  7. 53(10)”°  ^ — . For  actual  inputs 

to  Itract  the  edge  conditions  of  data  groups  148,  149,  207,  and  208  were 

used  to  simulate  conditions  downstream  of  the  shock  wave  of  group  150. 

This  led  to  an  Me  of  approximately  7.1,  a Tg  of  123.1  R,  and  a pg  of 
, lbfsec2 

1.1 7(10)”H j — . The  results  of  Itract  were  Included  with  those  of 

ft 4 

General  Electric  in  fig  10  with  tabulated  results  In  Table  XI  of  Appendix 
G.  Again,  the  results  showed  Itract  passing  through  the  field  of  laminar 
data  points  and  settling  high  In  the  fully  turbulent  region.  In  the 
fully  laminar  region  Itract  was  within  2.6  percent  of  the  data,  and  with 
the  exception  of  one  point,  Itract  settled  within  9.3  percent  of  the 
data  In  the  fully  turbulent  region  for  the  last  1.5  ft  of  the  cone.  For 
the  second  test  a cold  wall  experiment  was  considered,  data  group  1. 

The  same  free  stream  and  edge  conditions  existed,  and  only  the  Tw  was 
changed.  The  wall  was  cooled  from  1060  R to  580  R and  the  experiment 
was  repeated.  The  results  of  this  comparison  were  Included  In  fig  11 
with  a tabular  summary  In  Table  XII  of  Appendix  G.  Near  Identical 
results  were  noted  among  the  three  theoretical  codes:  Itract,  Nsbl,  and 
Vlzaad.  Nsbl  and  Vlzaad  were  codes  used  by  General  Electric  to  check 
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Martel luccl,  Laganelll,  and  Hahn  Data 


Station  (In.) 

Fig.  10.  Comparison  with  a Hot  Wall  Reference  Test  Case, 
Heat  Transfer  Group  150  (Ref  15) 
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their  experimental  results.  In  this  test  Itract  was  within  8.9  percent 
of  the  data  over  the  laminar  region  and  maintained  a consistent  13  to  15 
percent  high  prediction  over  the  entire  turbulent  region.  Consistent 
with  the  results  of  data  132,  these  last  two  test  cases  were  predicted 
best  using  eddy  model  one.  Having  noted  the  consistent  trend  set  in  these 
three  heat  transfer  cases,  attention  was  directed  back  to  profile  data 
groups  74  through  79. 

Having  used  the  output  of  these  groups  for  the  Investigation  of 
data  group  132,  the  profile  data  group  74  was  again  used  by  Itract  to 
predict  the  profile  shape  of  velocity  and  temperature  versus  for 
station  s equal  to  .466.  The  results  were  included  in  fig  12.  Due  to 
the  questionable  data  points  for  ^ less  than  .4  no  percentage  error  was 
Included. 

These  results  represented  the  best  predictions  for  heat  transfer 
obtained  during  the  study  of  the  cases  for  no  mass  transfer.  As  with  the 
flat  plate  study,  numerous  combinations  of  grid  size,  Prt,  and  eddy 
models  were  attempted  In  order  to  minimize  the  error  in  neglecting  higher 
order  finite  differencing  terms  and  best  describe  the  flow  behavior. 

Having  completed  the  cases  with  no  mass  transfer,  study  began  in  those 
cases  with  transfer. 

Beginning  with  the  greatest  blowing  rate  of  1.5(10)  , data  groups 

66,  68,  and  73  were  chosen  for  consideration.  It  was  found  that  Itract 
was  neither  able  to  predict  the  heat  transfer  of  data  group  66  nor  the 


nondlmenslonal  profiles  of  data  groups  68  and  73.  Various  grid  sizes 

v. 

were  attempted,  which  in  the  extreme  cases  Included  a An1  equal  to 
1.25(10)”*,  250  divisions  In  the  grid  along  the  streamwlse  direction,  and 


150  divisions  In  the  grid  along  the  n direction.  The  ratio, 


-,  was 


46 


decreased  to  a value  of  1.05.  Even  with  the  finer  mesh  size  the  tempera- 
ture change  at  the  first  two  stations  at  which  mass  transfer  was  occurring 
was  so  great  that  the  numerical  scheme  failed  due  to  attempting  undefined 
arithmetic  operations  related  to  these  temperature  differences.  One  step 
prior  to  failure,  the  coefficient  of  skin  friction  and  heat  transfer  were 
seen  to  be  decreasing  rapidly.  This  was  indicative  of  a numerical  separ- 
ation of  the  boundary  layer  and  the  imminent  failure  of  the  computer  code. 
A smaller  transfer  rate  of  1.(10)  nominally  was  attempted  next. 

Data  group  60,  depicting  heat  transfer,  and  data  group  59,  depicting 
profile  data,  were  chosen  as  test  cases  for  investigating  a mass  transfer 
rate  of  1.(10)  . This  was  the  first  case  involving  mass  transfer  In 

which  Itract  was  able  to  complete  the  calculation  of  the  boundary  layer 
without  terminating  In  an  error  finish.  This  did  not  imply  the  accuracy 
of  the  predictions,  only  that  the  finite  differencing  scheme  was  able  to 
proceed  through  a complete  computation  of  the  grid  of  nodal  points. 

As  the  profile  data  group  59  was  the  only  field  data  associated  with 

data  group  60  for  heat  transfer,  the  Information  from  group  59  alone  was 

used  to  determine  the  Inputs  to  Itract.  For  computation  purposes  Itract 

was  provided  the  following  pseudo-infinity  conditions:  M was  approxi- 

5e  lbfsec2 

mately  equal  to  6.7,  Ta  was  112.4  R,  and  pa  was  1.26(10)  3 — . From 

e e ft4 

the  graphical  presentation  of  St^  versus  station  along  the  surface  of  the 

cone  an  initial  transition  point  was  chosen  to  be  over  two  ft  from  the 

tip  of  the  cone.  Also,  from  tabular  and  graphical  presentations,  the 

ratio,  T^T^  was  approximately  .57.  Related  to  the  blowing  rate,  the 

supposed  actual  rates  of  transfer,  Xm,  were  8.3(10)  from  9.5  In.  to 

22.  In.,  8.(10)"*  from  22.  In.  to  34.5  In.,  9.6(10)”*  from  34.5  In.  to 

47.  In.,  and  9.(10)  from  47.  In.  to  the  end  of  the  model.  This 


(pv)w 

disagreed  expectedly  with  the  figure  for  y-jy-  from  data  group  59  which 
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was  6.3(10)"’’.  Initially,  the  blowing  rates  for  x were  chosen  for 

oo 

testing. 

Initial  testing  with  the  aforementioned  inputs  led  to  a series  of 

error  finishes.  Itract  was  able  to  compute  for  the  first  3.5  ft  of  the 

cone  at  which  point  the  coefficient  of  friction  and  Stanton  numbers  had 
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decreased  rapidly  to  values  of  10  or  10  . At  this  point  Itract  simul- 
ated boundary  layer  separation  with  an  error  finish.  Again,  many 
combinations  of  grid  spacing  were  attempted.  The  transfer  rate  seemed 
clearly  too  great.  With  the  lack  of  clarity  of  a transition  point,  an 
attempt  was  made  to  run  the  program  assuming  turbulent  conditions  from 
the  tip  of  the  cone.  With  this  one  change,  Itract  was  then  able  to 
successfully  solve  the  boundary  layer  problem,  but  with  two  conditions 
at  input  still  in  question.  First,  further  scrutiny  of  the  heat  trans- 
fer curve  showed  justification  for  choosing  a transition  at  1.5  ft  from 
the  tip  of  the  cone.  Then,  to  be  consistent  with  the  newly  defined 
pseudo- Infinity  conditions  downstream  of  the  shock  wave,  the  proper  mass 

(pv)w  (pv)w 

transfer  rate  was  thought  to  be  ^eu  • From  the 

transfer  reading  of  data  group  59  a scaling  factor  was  used  to  adjust 
the  blowing  rates  from  8.3(10)"*,  8.(10)”*,  9.6(10)"*,  and  9.(10)"*  to 
5.4(10)"*,  5.3(10)"*,  6.3(10)"*,  and  5.9(10)"*  for  the  four  sections  of 
the  cone  previously  mentioned.  With  these  adjustments,  Itract  was  again 
run  for  the  final  test  of  data  groups  59  and  60.  The  results  of  the 
heat  transfer  study  were  Included  In  Table  XIII  of  Appendix  G with  a 
graphical  depiction  In  fig  13.  In  the  turbulent  region  Itract  over- 
predicted the  experimental  heat  transfer  data  by  about  70  percent  with 
a 30  percent  average  In  the  laminar  region.  In  the  profile  results  of 
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Fig.  13.  A Heat  Transfer  Comparison  for  Data  Group  60  (Ref  14) 
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fig  14  there  were  Identical  temperature  predictions  near  the  wall  with 
a disparity  greater  than  46  percent  near  the  center  depth  of  the  boundary 
layer.  Concen^^*b^_"‘jlQ34j^jv£f*le..  the  General  Electric  data 
depicted  a near  separated  condition  at  station  s equal  to  .658,  and  the 
disparity  greatest  ?»ear  the  middle 

0 • ft* . 

of  the  boundary  layer  thickness  with  a 58  percent  error.  Observing  the 
near  'll'  ifes-STNe  low  heat 

transfer  data  results.  This  concluded  the  final  Investigation  of  data 
groups  59- and  60.  One  final  case  with  a nominal  mass  transfer  of  5.(10) 
was  then  selected. 

From  experimental  results,  data  group  203  was  chosen  to  study  heat 
transfer,  and  data  groups  200,  201,  and  202  were  chosen  to  study  the 
profile  characteristics  of  the  boundary  layer  for  this  lowest  mass  trans- 
fer case.  Free  stream  conditions  included  an  equal  to  8.0,  a T equal 

_5  Ibf-sec2 

to  98.1  R,  and  a p of  7 .48(10)"° » — . From  groups  200,  201  , and 

ft* 

202,  the  Inputs  to  Itract  for  the  study  of  group  203  and  the  heat  trans- 
fer consisted  of  an  Ma  equal  to^approximately  7.1,  a equal  to  120.6 

e «• — 2 e (pv) 


R,  and  a of  1 .18(10) 
v-4 


4 lbfSec 

' ~7F 


VTc 


was  .78  and  a constant 


(pul 


w 


equal  to  3.1(10)  was  used  as  the  transfer  rates  computed  at  the 
three  profile  data  stations  were  nearly  equal.  The  results  of  the  com- 
parison between  Itract  and  the  experimental  data  of  group  203  were 
sunmarlzed  In  Table  XIV  with  a graphical  presentation  in  fig  15.  There 
were  no  laminar  data  points  with  which  to  compare,  but  in  the  turbulent 
zone  Itract  underpredicted  the  heat  transfer  by  a 30  to  50  percent 
margin.  Noting  the  sensitivity  of  the  code  to  even  small  changes  In 
mass  transfer  rates,  data  group  203  was  retested  for  possible  actual 
mass  transfer  rates  of  1.(10)"*  and  2.(10)"*.  The  numerical  predictions 
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Fig.  14b.  A Profile  Comparison  for  Data  Group  59  (Ref  14) 
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were  shown  to  pass  through  the  region  of  turbulent  data  points,  also 
shown  In  fig  15.  The  study  of  data  group  203  provided  theclosest 
results  of  itract  for  the  investigations  that  included  mass  transfer, 
and  the  corresponding  profile  results  of  data  group  201  were,  likewise, 
the  best.  A comparison  of  Itract  with  the  profile  data  of  station  s 
equal  to  .646,  data  group  201,  was  included  in  fig  16.  Near  a ^ of  .1 
the  temperature  profile  was  33  percent  in  error  with  a 20  percent  error 
in  the  velocity  profile  for  a similar  boundary  layer  depth.  Both 
error  figures  represented  the  extremes  in  error  between  the  numerical 
results  and  the  experimental  data. 

With  this  test,  the  investigation  of  the  cone,  both  with  and  with- 
out mass  transfer  had  been  completed.  A summation  of  the  investigations 
of  the  cone,  as  well  as  the  plate,  followed. 
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Fig.  15.  A Heat  Transfer  Comparison  for  Data  Group  203  (Ref  15) 


V.  Summation 


Originally,  the  computer  program,  Itract,  incorporated  a boundary 
condition  of  zero  mass  transfer  at  the  surface  In  calculating  the 
boundary  layer.  With  the  program  modified  to  accept  the  condition  of 
mass  transfer  at  the  surface,  boundary  layer  flows  perturbed  by  this 
mass  transfer  could  be  solved.  The  purpose  of  this  study,  then,  was  to 
modify  the  basic  code  and  verify  this  modification  through  comparison  of 
the  numerical  results  with  analytical  expressions  and  with  published 
experimental  data.  Data  was  chosen  from  experiments  on  both  a flat  sur- 
face and  an  axisymmetric  cone. 

From  the  study  of  flow  over  a flat  plate  four  results  were  out- 
standing. First,  the  grid  size  was  of  fundamental  importance  in  solving 
the  problem.  A finer  mesh  of  nodal  points  yielded  better  results  to  a 
point  where  the  effects  of  truncating  higher  order  terms  in  the  finite 
difference  expressions  became  insignificant.  Second,  the  cases  investi- 
gated with  suction  were  clearly  more  stable  in  computation.  Further, 
these  cases  were  more  accurate  predictors  of  the  experimental  results  to 
the  extreme  of  the  asymptotic  suction  limit.  Third,  the  results  for 
the  blowing  cases  were  less  accurate,  and  the  error  did  not  show  regular 
trends  insofar  as  a fixed  error  amount  or  a fixed  percentage  error.  The 
heat  transfer  predictions  were  low.  Fourth,  for  the  case  of  blowing, 
the  best  results  were  obtained  by  using  eddy  model  zero.  However,  for 
the  cases  with  suction,  eddy  model  one  provided  the  best  results. 
Overall,  the  modified  code  was  verified  for  flow  at  an  M,,,  much  less  than 
one  over  a flat  plate.  For  both  laminar  and  turbulent  flow,  the  code 
was  proven  to  be  accurate  for  the  case  of  the  blowing  parameter  to  a 
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strength  of  1.(10)  . For  the  suction  case,  the  code  was  accurate  to 

the  suction  asymptotic  limit. 

From  the  study  of  flow  over  an  axi symmetric  cone  four  results 
were  noted  as  outstanding.  First,  the  grid  size,  again,  remained  an 
important  factor  in  the  success  of  the  numerical  predictions.  The 
finer  lattice  of  nodes  yielded  better  and  better  results.  Second,  the 
case  of  suction  was  not  studied  but  for  the  case  of  blowing,  the  pre- 
dictions became  erratic  as  the  blowing  parameter  was  increased.  The 
resulting  errors  did  not  show  a systematic  trend.  Third,  the  best 
results  for  the  cases  of  positive  mass  transfer  occurred  when  eddy 
model  one  was  used,  unlike  the  results  of  the  flat  plate  study.  Fourth, 
the  results  of  these  blowing  cases  were  shown  in  fig  15  to  be  extremely 
sensitive  to  the  blowing  parameter,  and  the  precision  with  which  the 
blowing  rate  was  measured  would  have  to  be  considered  in  completely 
evaluating  the  validity  of  the  modified  code.  Overall,  the  modified 
code  provided  reasonably  predictive  results  in  the  case  of  laminar  and 
turbulent  hypersonic  flow  over  a slender  cone.  Specifically,  for  a 
mach  number  of  eight  the  code  provided  reasonable  results  for  mass 
transfer  rates,  defined  as  (pv)w/(pu)e>  up  to  3.1 (10)"4.  To  verify  the 
code  within  an  acceptlble  limit,  the  precision  of  the  measurement  of 
the  blowing  rate  would  have  to  be  quantified.  Assuming  a measurement 
error  between  1.(10)‘4  and  2.(10)"4  was  possible,  the  code  was  verified 
for  turbulent,  hypersonic  flow  over  the  cone  for  mass  transfer  rates  up 
to  a strength  of  3.1(10)“4. 

With  the  limits  of  the  code  specified  for  the  particular  cases 
studied,  factors  that  contributed  to  the  obvious  limits  of  the  code  for 
the  positive  mass  transfer  case  Included  the  following:  First,  at  the 
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Initiation  of  blowing,  sharp  temperature  gradients  in  the  streamwise 
direction  resulted  in  numerical  problems  for  the  code.  Second,  with 
this  temperature  change  in  the  streamwise  direction  normally  considered 
insignificant  as  a boundary  layer  assumption,  the  effect  of  blowing  may 
have  violated  a basic  proposition  in  derivation  of  the  boundary  layer 
equations.  Third,  if  the  flow  were  separating  from  the  solid  boundary, 
as  it  seemed  to  do  in  some  of  the  velocity  profiles,  another  basic 
proposition  of  boundary  layer  theory  was  violated,  and  the  imminent 
arithmetic  mode  failure  of  the  code  was  to  be  expected.  The  success  of 
this  code  ultimately  depended  on  the  condition  that  the  classical 
boundary  layer  assumptions  were  not  violated.  Finally,  in  at  least  the 
study  of  the  conical  flow  It  has  been  found  from  previous  study  that 
though  It  was  valid  to  use  experimental  data  to  describe  the  flow  envir- 
onment downstream  of  the  oblique  shock  wave,  this  could  have  misrepresented 
the  needed  Inputs  of  this  code.  Further,  It  has  been  found  that  the  near 
adiabatic  condition  of  a wall  has  been  a most  difficult  problem  for  a 
finite  difference  scheme  to  compute  accurately,  more  so  than  in  the  cold 
wall  case  as  was  shown  in  the  favorable  results  of  fig  11  (Ref  19). 
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Appendix  A 


A Background  and  Derivation  of 
Some  Key  Expressions  Used  in 
the  Analytical  Solution 

The  differential  equations  which  described  two-dimensional  laminar 
boundary  layer  flow  in  a cartesian  coordinate  system  were 
Continuity 

,38, 


Momentum 


3u  . 3u  . „ 3u 

0 it  pu  lx  + py  5^ 


3y  3x  3 y 
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P 3t  + PU  33T+  pV  3y 
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Energy 
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From  Reynolds  the  following  definitions  were  made  to  describe  a turbu- 
lent boundary  layer: 


u+u',p=p+p',t 


yx 


T + T 

yx  yx 


(42) 


pu  = pu  + (pu)\  pv  = pv  + (pv)\  px  = px  + px 


1 fT+T/2 

T u dt 

' 'T-T/2 


(43) 


where  bars  Indicated  mean  values  and  the  primes  designated  Instantaneous 
fluctuations.  Finally,  the  definition  of  time  averaging  was  necessary 
and  was  explained  by  the  following  example: 

rT+T/2 
-T/2 

where  T was  used  in  this  example  to  represent  time,  not  temperature. 

With  these  basic  definitions  and  assuming  steady  state  conditions,  the 
laminar  equations  could  be  transformed  Into  descriptions  of  turbulent 
boundary  layer  flow. 

To  ultimately  reach  the  form  of  the  equations  listed  In  Eqs  (8), 

(9),  and  (10),  the  steps  were  Included  for  the  simplest  case,  continuity. 
Time  averaging  and  substituting  from  the  above  definitions  yielded: 


(pu  + pru‘)  + ^7  (pv  + p'v‘)  ■ 0 
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(44) 
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It  has  been  accepted  that  the  p'u'  term  was  strongly  uncorrelated,  and 
this  term  was  eliminated  from  further  consideration.  Then,  following 
two  coordinate  transformations  the  final  form  of  Eq  (8)  was  reached. 

An  assumption  of  this  study  was  that  flow  could  be  considered  two- 
dimensional.  Further,  a body  oriented  axis  system  was  employed  for  both 
the  flat  plate  and  axisymmetric  cone.  Finally,  a cylindrical  coordinate 
frame  was  chosen  to  describe  both  of  the  flows.  Performing  the  cylin- 
drical transformation.it  was  found,  first,  that  in  cartesian  coordinates 

3(PU)  . 4 ♦ *TT 

m 1 rv  . _ 


3x 


ay 


= o 


(45) 


having  dropped  the  time  averaging  symbol  from  the  mean  quantities.  Then, 
by  defining 


pu  as 


pu,  0,  p 


fv+^l 

P 


and  employing  the  definition  of  the  divergence  of  pu^  or  V-pj^,  continuity 
in  a cylindrical  frame  was  shown  to  be 


3(rpu) 

3x 


* 8(4 


♦ 23 


3y 


= 0 


(46) 


By  Including  an  exponent  with  the  r term  to  yield  r^,  it  was  noted  that 
by  setting  j equal  to  zero  or  one  would  yield  the  expressions  for  contin- 
uity related  to  the  flat  plate  and  to  the  cone,  respectively.  Then 
having  demonstrated  a transformation  to  cylindrical  coordinates,  it  was 
reassuring  to  show  also  that  a body  oriented  axis  system  x',  y'  could  be 
used  In  the  case  of  the  conical  flow  as  an  x,y  system  had  been  used  for 
the  flat  plate.  Figure  17  was  Included  as  a pictorial  description  of  this 
situation,  with  the  prime  symbols  serving  here  only  to  differentiate 
direction. 
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Fig.  17.  Showing  the  Equivalence 
of  Expressions  in  Rotated  Coordinates 


First,  it  was  recognized  that 


y = x'  sin<j>c  + y'  cos«j>c 


x 3 x'  cos4>c  - y'  sin<|>c 


y'  3 y cos<t>c  - x sin4>c 


x'  3 y sin*  + x cos<J> 
vy  - ux,  sin<j»c  + v , cos*c 

ux  = ux'  cos*c  " vy'  si"*c 


(47) 


Then,  from  fig  17,  it  was  true  that 


3(rpu  ) 3(rpv  ) 
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(48) 


If  F were  equal  to  (rpu  ) and  G were  equal  to  (rpv  ),  It  was  demonstrated 

A Jr 

that  the  chain  rule  could  be  used  to  ultimately  produce  expressions  for 
3F  aF 

and  such  that  the  following  equality  was  true: 

3(rpuJ  a(rpv  ) 3 (rpu  i ) 3(rpv  ,) 

-Tr*-*-V ST-  + ‘ 0 («> 

Thus,  through  two  transformations  the  immerging  expression  for  continuity 


was 
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which  matched  Eq  (8). 

In  the  same  manner,  but  with  increased  complexity  of  expression, 

the  equations  of  momentum  were  written  as  follows  using  the  equation  of 

continuity:  ? 

3(pu  ) 3 (puv ) ^ 3px  3fyx 

3x  + 3y  3x” + 3y 


3(puv)  3(pv2)  3p  3t 

+ = -*-£■  + — 
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(51) 


Employing  the  equation  of  continuity,  and  with  the  substitutions  of  Eq 
(42),  it  was  noted  in  the  final  form  that  7 was  much  less  than  7 and 
that  Eq  (51-2)  became  a negligible  expression.  Eq  (51-1)  was  dominant 
by  an  order  of  magnitude  analysis,  and  after  dropping  the  bar  symbol  over 
mean  quantities,  reduced  in  the  steady  state  case  to 


Then,  using  Eqs  (41)  and  (42),  discarding  negligible  terms,  and  trans- 
forming to  the  cylindrical  coordinates,  Eq  (52)  reduced  to 


pu 


M’ 


+EHI  ly.  = _ + i 


3y 


dx 


ri  3y 


(rj(p|^+  piTF-))  (53) 


which  was  the  momentum  equation  In  Eq  (9). 

Finally,  and  with  still  greater  complexity,  the  rules  of  substitu- 
tion of  Eqs  (41)  and  (42)  along  with  the  Idea  of  time  averaging  and 
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coordinate  transformation  could  have  been  employed  with  Eq  (40).  Then, 
following  steps  similar  to  those  of  Van  Driest,  the  energy  equation  could 
also  have  been  simplified  to  the  form  shown  in  Eq  (10)  (Ref  9:145-150). 
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TF(  J?0  A.  “n,  i)|  wpi**{6. ’"?H 

TF(  J?0?.*IF,  0)  N^UKS  ,*022) 

C 

P TNOIJ’'  T'lTTUL  OOCCIL; 

r 

1»  MS’-lOTr-* 

C INI'HinC  TM£  ST9;SMWTsr  LOCATION 

S«SI 

"S?  * C5  1 * "*S 

CX2CS*  ''XtnS*CXr>S*1. 

SE°C=1 . 

T,  INITI*  LT7I  HI  STRUMWTSF  LOCATION 

y(1>=o  ,n 

CC  291  IL*?,290 
IYsXXX  (LL*?)  *nyW 
211  Y(LL)sY(L.-H*OY 

CO  ’ll  LL  * I,  210 

ClILL)  *C?<LL>  * C3  (LL)  *X*IN{LL>  =0. 

7P(  LL)  ^VI(LL)  = Y0(LL)  =-yjll) 

rP<  tl»  *C1(LL>  = FN(LL>  *”’  (LL)  = TN(LL)  =TO(LL  ) *;»(  L L ) *i  3 (LL ) =£N  (LL  > * 

1 £TO<LL*  o-TO(LL)  sITN(LL)  *1.9 
7 9 0 CONTINIJ" 

CO  791  J * 1,  290 
CO  7flt  T ■ 1,  3 

791  Al<  j,T  )*!?(  J,T  >*A3(  J,  I)  «91  (J,I)  *32  ( J,I)*A3(J,I)  *C1<  J.I) 

1 «C2(J,I)*03(J,I)*9. 
oorcsO  • X'<I*lF» • 2 
Torp  > (1  - 1 , ) * YNlNr ** 2 

INIYTALI’i  COUNTERS 

IOOMN^-art 
To*  ICC  rc 
TG*  l 
IP*1 
n»Cw*  0 
TTCNT1  * 1 
TIN*0 

• ••  V»!m  FI9ST-0^C2P  TPIOIAGONAL  HATRIX  SOLJTycn  ••• 

CO  115  m*n5T1»y,ICn01 
TF(*t.rc,Nsr»oTi  vp.nst* >* 

TFIC.EC.  TSNCl)  ‘*«i* 

IFM.yc,  (n/NNP)»m«;S)  NO** 

5*S*CS* 

9X295  * 9X135 


093991 

103«11 

999920 

991930 

901340 

990959 

190359 

010970 

009959 

999899 

990911 

903910 

900920 

000930 

000940 

910950 

100969 

100970 

100930 

000991 

991199 

101010 

101020 

001939 

091040 

001050 

091069 

001079 

091031 

131090 

091110 

001113 

901120 

001130 

90111.0 

001150 

001160 

001170 

001139 

001190 

901200 

001210 

991220 

001230 

001240 

001259 

901260 

901271 

901230 

101290 

001300 

001319 

001320 

001330 

001340 

001350 

001369 

001370 

001330 

001390 

001400 

001419 

001420 

001430 

001440 

001459 

001460 

001470 

001430 

001490 


68 


c 

c 

c 


XI  39 


IX?1; 


n 

C 

c 


64? 


6?6 
6 58 


C 

c 

c 


l 


c 

c 


c 

c 

c 


c 

6998 

c 

c 


G0'*:>'lT"  !.o;«U  t '10  on.^5U?i  GRADIENT 

-a.,L  JOf^u(<;,  yMT.|P,*.  ,^a',  1 , 0 °9G  l f r?  t>|F  , X HE) 

CqmoiiT  - LOCAL  60G"  o’O’P’.TIiS 

3C  - 

50  = TCOGl  /P’’.' 
te  = r'f'i'/p;f 
IJE  * 3'»T(2.*<T10  - Tt>> 

9E=G*  = <v  ((5.i,  i))  **r) 

TP=SU/  <T2TNF»T  »> 

IF  (OH^Ol)  642,675,642 
XKMEsT  "*  *0"EGA 
GOT -362  » 

xlr.i"=T  r**i.  5*<  1.  *138.  V < Ta*rpF?))  / (TI*iq«,6/t  TA  «t?ef) ) 

c cn  r m "■ 

"OHO'IT'"  LOCAL  XI  BNO  CTE9  LENGTHS 


0Y3S=?f*'IE*XNII  r 

IF(  J?3  ?.  NF.  0)  OXOS  = OX  0I»S*»2 

TF(H.E'\?>  3Xns  = 0X?9S=9X0S 

0X?=.  5 *9S2*  ( (1  .*0$2/9Sl>  »0X1  D?*0S1»3X0S/  <DSl*3S?)  -0S2* 0S2* OX 20S/ 

nsi*  nGt*0S2)> ) 

3EYN4E  rPF*U2‘S/XNUE 

’>~y'y'r  =9fy*  vis  i»iF*qsr  »»d  =: 

IF<h.  ?'•.’)  OX  1 =0  X2 
TF(M  x=0X0S*5t 

x=x*ox  2 

OO'OUTF  STE°  LENGTH  c UNCTION s 

vi*».»  nxi*2.»ox?>  /( o xi  * ox  2) 

inni"  .ft.  u »i  i ». 

X?r((OXl  *0X2>  /OX  11*2.  0 

Y3s  (OX  2*0X2/  (0X1  ♦< 0X1  *0X2)1 ) *2,  0 

yl=  (ox  i* ox?)  'oxi 

Y5=OX2/OX1 
TWTE  = T'H/TE 

COH°t)TF  AL»UA,  3ET4,  INC  LAH90A 


OUEOX=  -=■>/(  =>F»  U2»9X0S  ) 
X4L=UI  *""/T“ 
x BE  a 2 . 0*x*OUEOX/ur 


L'NGTh  =IEOGE 


ASSIGN  THE  HSTOIX  ELF“-*‘ItS  fQ®  THE  FINITE  OIFFE^SIE  EGUaTIONS 


001500 

001510 

001520 

001530 

001540 

001550 

001560 

001570 

001590 

001590 

001600 

001610 

001620 

001630 

001640 

001650 

001660 

001670 

001680 

001690 

001700 

001710 

001720 

001730 

001740 

001750 

001760 

001770 

001790 

001790 

001800 

001910 

001820 

001530 

001540 

001950 

001360 

001870 

001380 

001890 

001900 

001910 

001920 

001930. 

001940 

001950 

001950 

001970 

001980 

001990 

032000 


'“ALL  *L«1TX(  H,0X2,X  ,Tf.L,X9E,T9,  IO I F* , y 1 , Y 2,  Y 3 ,Y4 , Y5 , T HT-E  , I TCNT 1 , 00  2010 

u,a2,  13,91,97,9T,ri,C2,23)  002020 

ASSIGN  the  N»T°IX  £LE«ENTS  C9F  THE  FINITE  9IFFf?en;E  EOUATIONS  002930 

002040 

002050 

"•T2IX  riVEASIOM,  SOLXE  FOP  F,  THETEA  AND  X 002060 

002070 


**L  H ACTONS  ( r p,  to,  >/o, 31,32,  03,  41, 9 1, C 1, 42, 32,  C9,  A 3, 33,  C3,  3,  LENGTH  00208  0 


,210)  002090 

002100 

«arst«  t*wf»sion,  sol r f«jo  f,  thjiea  ant  v 002110 

002120 

•T*«»TJ  t 002130 

002140 

Mrt*0*H***l  t (OA(*(l  002150 

/ if  « l.*t  ./x«<|l  -vPf  IFOGE-t)  •(  t. -t.  7XXO  *XX<-  002160 

+ • M ••«*/<  t,  •»./*«*»)  /BY  00217  0 

».  902160 

002190 


’r'i5T'l=M,<Vl 
~ V=  3V'J  **  w * * ( .J  -1  ) .IV 
CP(  » «|  xT  •>(  I')  a 1.  1 
*,5  »C>(  ’T)  = /■>(:- I ji-  1)  »>/<  *'V 

C ’M?TT1  T”M  ir  PIUTI_1P  '?U-ITI  VJC 

TFr<.r''. ’)  r.o  ro  "'jj" 

*o  *n  2919 

• not)  «TQ  «1io  * s 1 , < 0 M 

'/0<  !>  = ')M(T)  =V°  (I) 

-0(  XI:  r*l(T)  = eo  d ) 

•<n i<a  tc< I > * ”►» c I » sT®  (I) 

r rntr:& *-  in  or  sr"in~  solutions 

101?  *0=t"'f.<- 

c 

c IJ  SN'O  *M-TS  PROFILES  I'SpaTIONS 
-aj?=(  e-3  r>>  -cp  (i  )>/o»n 
TF(tT*  M*t.  CO.  •>  ) TaiJl=  11.  *TAU2 
3T1  ? = T "’JI/TS  J?-1. 

•’auis’'  ai » •» 

tFdT1**'*!  ,Lr.  100)  39  TO  70  05 
viot  T~ { S ’TO!)  ** 

08LL  rxt’ 

79 1e  IF<A1S  P*i2>.GT.  "290?)  OO  TO  5959 
o u a *n  ’■wcra  ppofilps  :*f=ationp 
c 
c 

c 00MOIITC  OLT,  O0T(0*LTl  STAP)  AN9  34  T (THETA) 

5-  OP=T3(t) 

TPT*1. 

2LT  *3L  PT*PL  4T  = 0. 

TWdl  =0. 

00  5’  *'*?,X0*I 
9Y*QYV»XY<*»  (N-2) 

* *t  p ( *n 

ToisTPTt.,?»0T*  (CO+3) 

on=c 

y)ri<N)  =t->t»so°t<?.*x> 

If<  J?o  a. ‘15. 0)  x*iu(‘i)  * xmn  m / s 

tL0Ts9LTp*(?.  - CP  (5)  /*P('1>  -fo(M-I)  XT°<M-l) ) • <X‘IN(‘I)  -XNM  <N-1) ) /*. 
3L4**3  l'JT+  (FP(  N)  Ml.  -FP(M>  ) / TP(  *1)  +FP  <N-1>»  Cl.  -FO(N-l)  ) /TP(M-1 ) ) 

1 MXNN  <>l)-XN‘t<  *1-1)  ) /?  . 

TFIUT.OT.O.J  00  TO  57 

tF(FP(  t.)  .r.F.o.  gq5)  3UT=X»1»|(N)-(  FP  (»l)  -.  995)  • CXNN  (5)  - XNN  ( N-l ) ) 

1 /(r»(  »') -CP(N-I)  ) 

5*  "C^TTMI'T 

aLT  = 3.  T»ro5 
OLOTsO  1.0  ”*  • »5 
3L4T=3C.4T*rPS 

C IfrtoilT*  3LT,  9 OT  (OSU T 1 MAP)  ANQ  34T  (T mft  A) 

C 

C “COPtJTC  TMF  "00Y  YXSC  OSI^Y  C OECFICIp'IT 
TFCS.l  F.  OT’X)  GO  TO  5* 

rail  PFVOTP  (K0N,TP,X,'»2F,XMU£,X9£,S,ITCNTt,TP*) 

C OO-P'ITF  rHE  EOOY  YISOO'I'Y  C0EFFi;ir,,T 

c 

51  TT3NTI  *1 

C 

C 1SS55PF9T  or  GPIO  P0‘<ITS  1*1  ETA 

C 

TF(T»in'-M)  71,  71,  732 

71  QOMTTNU* 

IFM  - 71)  702 , 737,  77 

72  IF(A35  <r’D<r*0Ge-13)-co(:-0G£-l»)>-0,  9001)  70,73,71. 

73  TFU95  (TPIXSOGE-15)  - TP  ItCOGE-  16)  ) - .0001)  732,  732,  74 

74  I^GE*  ms-  *t 
IQ*  TO*  1 

9Y«OVM»yxK,»  (IFOGE-7) 

YfTPor.T)  * Y ( I FOG?*l)  *•  OY 
732  10  * 10-1 

C ASSiSl^lT  OF  GPIO  OQYITS  IN  ETA 

4 


092200 

992219 

902229 

002230 

902241 

902259 

002260 

902270 

002299 

902290 

902300 

002310 

002329 

002330 

002X1,0 

902350 

002369 

002370 

002360 

002390 

0021.00 

092510 

0025*1 

002930 

002540 

902450 

002460 

002470 

002490 

002490 

102500 

002510 

002520 

002530 

992540 

002550 

002560 

992570 

002560 

192330 

002601 

002619 

002621 

002630 

002640 

092650 

002660 

002670 

092660 

002690 

002*00 

002710 

002720 

002730 

002740 

002750 

002760 

002770 

002*50 

002790 

002600 

002510 

002920 

002630 

002940 

002659 

002960 

002970 

002660 

002690 


70 


c 

002390 

c 

002310 

r’  \ 

c 

^*veif  rr  S T:>  i *P  •''/if  TPMSf:»  5%|?  TjTO"~  STATTO'I 

032920 

*au.  rr-t^  c j*  »*#«»  *■**»- #T„}T,*iL*TOLrf  9r,ioai'tf?*: /rx-t 

002370 

1 X*  T « * p' 

002340 

c 

~C'f’'ITr  nAll  ST®"SS  IT  h?aT  ™a*iS"?  AN">  OUTPUT  STATION 

302350 

c 

002360 

r 

002970 

r 

SHIP”  o=1C:l*S  9A0K  ON*  XI  STATION 

002990 

r, 

002930 

NN  - I " * 5 

003000 

■>0  114  N«1,*IN 

003010 

FN(N)  = rOC|) 

003020 

cO(  N)  . rc  <*M 

003030 

TWM)*TMN) 

003040 

t o(  *n  i *ai  >n 

003050 

n»!(  *n  * wen 

003050 

V9<N).  tfP(N) 

003070 

*T-|  Cl)  =““9(N) 

003030 

“T9  (Ml  =-  ~=  CM) 

003090 

:M<N)  = <^cn 

003100 

cC<N)s“aCI) 

003110 

it* 

scntint 

003120 

5Xl*0<  » 

003130 

1S1=7S? 

003140 

TfCIH-'HSItOI)  114,113,114 

003150 

tl* 

OS?*’.  o*7Si 

003160 

IG  * ir.*t 

003170 

TNOCH  = 1 

003190 

TP  («.  aa.IENni  ) GO  TO  ?!’ 

003190 

GO  TO  ttl 

003200 

114 

OS®  = a3  ( 

00  3210 

TNI  OH  = 1 

003220 

TP  (H.  ra  • I E NO  1 ) GO  TO  237 

003230 

I f.. 

GO  TO  111 

003240 

l • 

2V 

t IN  * i 

_ . _ 

00  3750 

til 

"ILL  9®N5hS  (I00JN,  1°, 15,10, -START  ,IIN,H,5,V,3LT) 

003260 

115 

0 ON  T T N 1<e 

003270 

S'")* 

003290 

“NO 

003290 

SU990(J““NF  p»ESSM(S,XH,  5,0,1  onX  ,T,  YD 

003300 

rOHHOS  /O’OSTl/  0P(24)  ,<3(?4)  , ns  (24)  , losr$ 

003310 

111 

C0°H4T  f 5X, 'HAWING. ... CAtC'JL#TIO‘l  IS  O'JTSIOi  Oc  US 

®9ESC9I3E1 

99003320 

IPS?'!®!  0 OT  S , S IS  LESS  THAN  XP(1)»> 

003330 

?ni) 

eNHA“  ( “X  , • HA  WING.  . , . C ALCULA*  TON  TS  nuTSIOS  OP  t-|£ 

®9es:9I  on 

99003340 

l'SS'JO"  OITA,  s IS  50-AT-o  THAN  XOJtNO)*) 

003350 

390 

en’HAT  n X,SS15  .0) 

003360 

10.0 

003370 

T0H1  = I no'S-l 

003390 

IF{  T®®  rS,  PO,  0)  GO  TO  40 

003390 

00  ?9  t=i,:»»ps 

003400 

TP(  S.L  X=><  1 ) > H9ITE(6,190) 

003410 

IP<S.5'r.XOtTSt)PC))  H4ITS(6,200) 

003420 

TP(S.Lr.  T°<1>  ) I°=l 

003430 

TFCI9,  ‘IP.  9)  50  TO  39 

003440 

IF(S.;P.XO(IPHi)  ) jo.  I°PES 

003450 

IP(  r®.  Np.o)  50  TO  30 

003460 

1 1 

IFt  (S.  GP.X»(I)  ).  ANO.  (S.  Lt.X®  (!♦  till  19*1 

003470 

Ei 

TFII®,  po.O)  50  TO  20 

003490 

c 

KEEKTNG  tms  9pST  PIT 

003490 

TS*  <S-X»(I>  ) /(  XO  (1*1)  -<9(I>) 

003500 

TF<9S.  OT.O.  *)  I».IM 

003510 

c 

SEEKING  Op  ST  PIT 

003520 

IFII®.  NP.9)  GO  TO  30 

003530 

20 

OONTTN  ,c 

003540 

30 

IFCI».  GT.jphi)  IR.IPHl 

003550 

ir»*:»*i 

003560 

• 0H,I»-1 

003570 

_ 

Tp(  TB.  r®  • 1)  I®-*T?»2 

003590 

t 

C 

Th?  0U9IC  SOLIN*’  COS  PFI TINT S 

003590 

71 

c 


n 


40 


sc 


no 

200 

300 


c 

10 

r 


c 


20 


XI*  (yo  (•■’O)  *yo  (IS«)  -2,i*yo(T5))*(xp{IoP)-XO(I9*)>  003510 

V>-  (V3  f-aM)  - X3  (I -on  * ' X'  ( V“l  -X=  (T’l  ) 00  361  0 

X3=  (X3  <T’>'M  (T3)  )*(  X=  ( -XP  (I3>  > 00  362  0 

Y6a  to(  .»3(  I?3)  003630 

X5*X  = ( T3  ') -xo  ( I3.)  003663 

X6=X3(  ’•’»)•<’(  I?)  003550 

nrrs=<  r*  X6»X4  303660 

0 2=  (33  (I?)  *X1  ♦ 0°  f I53)  *»?-OP(  :®M)  «V3)  79  ETC  033670 

0,7=  (0=  (T?l  *X4-0O(l-»o)  »yc*no(  to«)»x5J /<5JTO  107681 

co“°mt ~ ~u~  r'.joir:  5=»l i**"  :n=  pft  cents  ooi69o 

nyos<;.  yo(  J = ) 003700 

lyoo-'j  vo*»»  003710 

->X3C=oyo '70.  00  3720 

0P9X1*  r>°  ( Is)  003730 

0="D(Io)  003760 

00  10  T=l,20  033750 

X=I*7X==  003760 

X2=  X*X  00  5770 

70'<X'»: ( 7s)  *3  2‘ x*7  3*  X’  003730 

3=p*0.  e*  (93nXl  *9P3X2)  *7x=F  00  3790 

3O0Xt="3TX?  003300 

TP3X=3°(I3J  *“2*9XO*C3*3xP7  00  3310 

T*s»**<  fC-1.  0) /G>  003320 

?.  "•{  ( 2.  0M5-  1.  0)  *XH»XM)/{?.  0»T>  -1.  0)  / Ol.  0)  1 00  333  0 

N9I7.=  ( ft,  300)  S ,o-.P=0X  YX  003840 

50  TO  c0  003950 

7=1.  0 003860 

OP0X=0  . 003370 

TrO»»(  (5-i.  0)  /G>  003880 

YM=S09"(7.1*  ( (2.  0»<S-  1.  0)  *XM»XM1/  (2, 0*T)  -1.0)  / (G-l.  3) ) 00  399  0 

9=TtJ9N  00  390  0 

= N3  003910 

SU9?0'J *TN=  SNT  '-*op(9T?  X,  3P«AX  f 5,  XNS0 ) 003920 

CON 40N /C39&TS/  CP(?4)  ,X°  (26)  ,0°  (24)  , 100=5  003930 

F0®M4T  flv,»  = I05T  30  ->s-s  eoiT  YIELDS  ADVERSE  PRESSURE  GRA0IENT  T3  003940 


10  3T"i°  =0’  C4LCUU4TI  0*1  *0  CONTINUE*) 

= 09X4  7 (1M0 , 11X  ,.3X5/U,  l'.X,2HC9,l  1X,5H9/9I-IF,  l’.x  , 4*00 OX) 

= C9XA7  (1 X . 6 ( 6 X ,'13.  9)  ) 

9°T  0L=  P0  9ax»l,  01 

roxoijrr  -ur  TRAILING  £03=  9=3X 

!PX1  = I 99*5-1 

TPN2  = I 997F-? 

1X1  =X=  (T3*U)  -XP(I=R  = 3 ) 

0X2  = X®  (7993) -XP(  1=9  = 3) 

0X1 3=9  xi  *0X1 
0X23=9  v?  *9X2 

np(  i3^c^)  = (Ca(  Tp *2)  •oxiT-CPt  Io**l)  *0X  32 -CP ( 1=’.* S)  * ( 3X12-0X22)  » / 
1 (9Xl,r,X9*(0Xl-9X2)) 

onxpuTF  th=  trailing  eige  ooox 
r*rtx=o 

CCXPUT*  '«=  L£  40 1*15  £355  3=9  X 
7X1  *X3  (9)-x  = (l  ) 

0X2aX=  (T)-XO(l  ) 

9X12=9 V1‘9X1 
9X?2*0X3»0X2 

9P(1>*  (C9(3)  *3Xl’-3=(  2)  • 9X22 -CP  (1)  • (3X12-3X2?)  ) n 9X1*  DX  2*  ( 3X1  - 
1 OK»l) 

IF(0°(  1)  .GT.OoxiX)  H9  IT  = (6,100) 

IF(3=(  1)  . IT,  3PMAX)  CALL  **TT 
CO*l“UTF  *M=  LE  40 IN'  £0*-=  3P0X 
30  20  T«’,IP*U 
TH1  = T-1 
IP1»I*  1 

9X1  «X=  (TNI)  -XO  (I) 

9X2«X9  (7=1 ) -X®  (I ) 

3X12*9  vnoxi 
9X22*0X7*9X2 

9P<  I)  = (*9(  IPi ) •9X17-;  pc  Mi)*  0X2  2-CP  (I)  *(9X12-9X77) ) /( 0X1  *0X2* 

1 (9X1- px-*)) 

IF(  (90  ft), 57.  0PT0L)  .A*r.  (X°(  I).  IE.  9TRX) ) INAX.I 


00  395  0 
003960 
003970 
003980 
003990 
006000 
004010 
004120 
006030 
006060 
004050 
004060 
004070 
OOhOOO 
004090 
004100 
004110 
004120 
004130 
006160 
006150 
004160 
004170 
004180 
004190 
006200 
004213 
006220 
004230 
004260 
004250 
004360 
004270 
004280 
004290 
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GO  T*>  *1 

C t«?  ro  utl  * *4  the  lEAGt-ig  EGG*  ""SION 

r ***  i x J «iv»i 
!*1:;  “1  »»1 
'xi  - to  r <‘iu -»»( tmaxi 
GYPsT3 

1X1  70X1*1X1 
'X7’=iy'*gx2 

l"P(  T^l)  *CP(IHD1)*0X  l? -CP  (I  '••4X  >*(0X  12-9X22)  OX  1*9X2*  <3X1- OX?) 
i •toms  Yi  rif>7 
GO  TO  19 

C S*OOthtnG  TMj  GO  UTS  JV|  THE  LEADING  EDGE  region 

«i  m<?:t-(  f-,  -"jot 

00  GO  I*t,IOEES 
°c»?.n »» '»( i)  - 1.  o)  /<  g*»«si> 
jn  hptt?(  «,  tg?)  xo( :)  ,»c <u , opt  i) 
eetupm 
f>r 

SL^EHOT'IE  crgriio  .XM'IE.X  , S , XEE , 1 ,8L9T  , 3L«T  , T , o9G  1 , 309  G 1 , 

1 P;V-if-, 

IC^NON  G,  30,  =ry,  XMI’I*.  PMEC-A,  00,  TW,  31.1),  Til,  0,3,  v IdO,  TE, 
1 BE,-5:  .'JE./rsrSP.SU.EP'.iS.DYW,  ci.EEEOE.TS.TA,  IE0GE.IEN01,  I*r  ACT, 

7.  PPT.X  w,ar?»,  XL  AM,  /l  pop*,  XI  NT-  =>,G£PO,  IC C “ > , T =9*1  <9>  1 El  (200)  , 

0 F>l  (00  "J  <?11)  ,rro  ( 30 0 » ,rTN  , E TP (210)  ,r9  ( ’GO)  ,**1  (210) , J23A, 

6 FP(’0  n)  ,TM(?1  0)  ,T->  (’  01 ) , XNN  <21  0)  , VI  (2 OH  , «(?1  0)  . X®(?0  3) , T=<200>  , 
5 01  <?o  01 ,03 (•>«,  o,  ,o3(?  oo) 

?11»  e02**AT  («41,11X  ,5HX/L  =,215.3) 

2101  'O®'**!  (?X,THXHE  = ,E1t..3,2X,7HP.E  = , E15.  9 , 2 X , t*Ooot NF  = , El 3 . 9, 

1 7X,7-(V0E  =. E15.3,2X,7mTW/TE  =,E15.3> 

Pin?  FORMAT  f?X,70L*  = ,E  lr.  3,2X  ,7U!>1*'  = ,E  1 =5 . ■ , ?X  .7H9L0  T = ,215.3, 

1 2X  . orvy  f - , 51  6.  3 , ? X,  7MBEY  OT  =,215.9) 

2003  p^NAT  (?X,7h:f-I0  3,  * I-1. 3, 2X  .T^GPE'IO  =,E15.  3,  ?X,7M5T.X0  = ,E15.8, 

1 ?X  , t*  PT  END  =,Elo.9,»X.7HPfYEX'=,E15.3) 

THTE=t  3(  d 

TF(OHEG1.£0.o.  ) GO  TO  >55 
TFtOMto&  .30.  1.)  GO  to  9551 

xlh  1 = i./(T3te**(i.  - omega) ) 

GO  ”0  »SG 

115  xlxi  = ( 1.  0»TP  ) » 50®  T ( T4TE)  / ( TWT  E*TE) 

GO  T0  1*5 
3551  XL-1  * 1. 

155  CONTl^ije- 

Yll  = <<  ?.  *XX<)  • (1.  *XXK  *-»X<*i’  ) *1  «f  XXO/  ( ( l.  »XX<)  •<  1.  *XX'<*XX<»»7)  ) 
Y12=<1  ,*XXK»  XX  «•••?>  /XXX* *J 
Y13*<1.*XXX*YX<*»7)/(XVY*»t»(1.*XX<)> 

Y1U  = 1./(XX<**3*(1.*XX<‘XXX**2)) 

TAU  = XL  “i»E£*XNU**,JE»'JE»  (-Y11»FP  (l)  *Y12*F3(2)  - Y13*P3(3)  ♦ Yll»*F?<5)  ) 

1 /(I**  *$1ET  <2.  »x  > ) 

05  * X L*U*EE*XNU2*UE*  TE*  ( Y 11  *TP  <1>  - Y l2*To<  2)  *Y  1 3*  T 3 ( J)  - Y14*Tp  ( 6)  > 

1 /(  1YH  *5 OPT  <2.  *X  > *°E  > 

TF<  1*0  S.  >15.  0>  TA  tl= T A 0 *c 
TF(  JPO  A.  'IE.  0 ) 05=05*5 
5TJ0  = 1, 

tFMO  .ME.  1.)  GTNO  = : °5*05  / (( 1.  - »0)*(TE  ♦ ,5»UE**?>) 

GTE  *10  = GT>n/(3E»JE) 

'FMO  * p.*EbS*  TA U 
CFE‘10  = 0**10/ < pf»UE*Uc) 

9FY0T*  p-yfxT’ILDT 
9FYMT*  P*YEXT*1LiT 

c 

C 5O.E0TT1M  0*  T*c  O'JTPIJT 
IF(H.N*.MP)  GO  TO  1811 
5tL ECT  T0*J  0*  TME  0*JT»'jr 

OUTPUT  STATION  OATA 

WRITEr5,2110)  5 

«PIT5(  ‘.'Oil)  XMF.P5G1, 0O9G1  , X9*t  T:rt 
MOTT* ( 5,  ’81  2)  0LT.9L5T.OL0T,  PET  IT , EE  YOT 
v«KT*«  5,’IO.T*  CFNO,  3*  E'lO , 5TM  0 ,S  TEMO,  EF  YE  XT 


00*300 

001.311 

011.321 

101*330 

001.31.0 

004350 

104350 

004370 

004330 

004391 

004400 

004410 

004420 

004430 

004440 

014450 

014460 

004471 

004430 

014491 

004510 

014510 

004520 

004531 

014540 

004550 

004560 

004570 

004580 

004590 

004600 

114610 

004620 

014630 

004541 

114650 

004660 

004571 

004690 

004690 

004700 

104710 

004720 

114730 

104741 

104750 

014760 

004770 

104790 

004791 

304800 

004910 

004920 

004830 

004940 

004850 

004860 

004370 

004860 

004890 

004900 

006910 

004920 

004930 

004940 

006950 

004960 

004970 

004991 

004999 
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I SET'PS 
"*0 

riai-ijrir  ri_MSTY(  * 

l ttomti,  Ai,A2, 

, s>,  j=v,  <« 

1 «•,'->*  .•I’.YTTIM’,’ J,  f 

2 P9T,  * vv,3T?<  , YLA«, /A 

3 r>j(’"  «)  , -3  r>11)  ,=T1( 

4 r=(’1  n>  ,'3(’ 

5 r»j  ,•»»  P)  ,1?<->1  1)  ,1T  (? 

?i**-*ii to-i  n (2 on ,=> , s 
1 «?(  21  ",  Tl  , Cl  ( 230,31  , 
COM'*0  “/’L-OAT  i/N'JHOS 
3 TIL*.  *L‘TM«T,  =V°AT 


,">Y?,y,  <AL,X3’,T=,IOI=e-,Yi,Y2,Y3,Y'.»,Y3,THTI  , i 

i’,n  ,32.3x,ci,C2,C3)  ' 

tv c t rt~r  A,  ’3,  TM,  =1",  T,n,  , j rr, 

«=.  •>$, oyw,  si,:<?90o,'c,ta,  T intact,  1 

■5^9  T,  XI  (T’9,5’90,  ICHS(S)  , t”N<9)  ,£Ot?nO>  f I 

2"") ,"TN(?1fl>  ,’-o<210>,r0(’l9>  ,rN ( 201 ) , J20A , 1 
r > ,XNN  (2n1),VU(293),V0(?3T>,Vo(?30J,T9{203), 
13 ) I 

?(?11,’),A3<?10,T),11<?00,3),32<201,3), 

C'(?<SP,  31,03(2"'), 3)  1 

T,X“OS<2fc>  ,=»H0V’AT(?(.)  , TOO  , T2  ;F,  XNU’ , 3 , 3S  , S t , 


*H*  T N ’ •?"  OfvjNOAlY  COM3  ITI  CM 
10  111 t T=  1 , 3 

t U(  1,T  1*12(1,1  )=S3(1,  X)=91  (1 ,1)  =32 (1,1)  = ’3(1,1)  =01  (1,1)  =02(1,1) 
1 *03(1, I>  = 0. 

£ 1 < Til  > = 1.0 
’2(1,1  >=1.0 
•on  1 > = ", 

02(1>  = TN” 

TF(S’t>«. ’0.1. 1 oo  TO  3012 
03(1,1  > = 1.0 
<C7Ws". 

1=  (T  OLO.LT.  1.AN0.S.  G’.SS.ANO.S.L’.  SO) 

6 CALL  T‘0*PLK(9<'V«,y,9’Y,  J204,Co=,9  = ,tjr) 

1=  (I  "LN«  GT • 0 • AN?. S.  Ci.SS.ANO.S.L’.S?) 

«i  CALL  ’•£  MOL M ( 0 CV  V,X  , * ’Y  , J?OA  , EPS*  , L?  ) 

03(1)  =10  VM 
30  TO  «113 
? XL=0X2 /( ’,  9»0Y  W) 

A 3(  1 ,1  ) = ’X2  *X  » Y1 

03(  1,1  > = -?.*  XL  M2.  MX  *)  / (t.*XX<) 

03(1,2)  = ’.  • XL*  (1.»XX<)/yxk 
*?(  1 ,3  ) = -2.  ’XL  /(  XY<*  ( 1.  »XXY)  ) 

0?(1>=". 


THE  TNM’»  tOGc  90UN0A  °Y  CONIITTON 

THE  FIFLO  ’OINTS  ’VALUATION 

I *IH1  = II  "3’-l 
00  SOI  4 •!=?,  NM1 
OV*XX<“  (N-l)  • ?YW 
0YN1*?  V/XXY 
XL»  OX?  /(  ?,  9 • OV  ) 

Y6*?./  (l.*OYH  /1Y> 

Y7«0Y/oyni 

Y«*2./  ((OYYl/O  V)  •(!.  » OY  “ 1 / OY  ) > 

Y9=  2 , / (t.*Ov/0YNi) 

Y 10  = 1 . -OY/OY-U 
SF®*1.  n 

TF(FO(M>  .LS.1.  ) S’°=1  • 

TF(  TTC  NTl  . GT.  1)  30  TO  ’030 
TF(  I?I  F=  .=0,  1)  GO  TO  ?501 

cmj  = Y(.»ro(‘J)-Y3*FN(N) 

TM1*Y(.*-")(N)-Y5*T'I(N) 

VH1  «Y4  •VO  ( N ) - Y 5*  VN  < *1 ) 

TF(  3=00. ’0.1.)  VU«XO(*n 

«H1  « (Y4»(’0(  M-l)  *10  (M)  +’0(M*1  ))“Y5*(FN(N-1)»EM(M>yEN(N*1I))/3. 

’Ty i *(Yfc*  (£TO(N-t)  »aT0(‘l)»£TO  (NYt)  )-Y5»(£T*((N-l)  >£TX(H)Y’TX(N»1 
l ll>/3. 

GO  TO  Till 
= H1  * F1(M) 

TNI  « TO(N) 

VN1  a UP  (Ml 

*Nl«(rO(N-l)  »£0(N)  *’0  ('!*1U/T, 
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1 


( 


’TUst  r*’{1-t)  **'0(9'  *-'T1(1*l ))  /?. 

GO  T">  'n’t 
7nnn  '•“l  = 

TMJ  i ”!>(11 
1*1  T:  '»’{*ll 

THt  * (*  *l-t ) *E  = (N1  *?’  fi»l)  )/  3. 
tt-i i t ~- ■» ( n- n *?Tnci)  /i. 

7101  I"  (11  'It  .:1.  0.)  -.1  *n  6 94 
nto**-:  gs  , pi  i.)  go  " *.6*- 1 
VL^l-l  • / { 7*11**  (1  .-01-  G1  ) ) 

YL’lls  (0  *734-1 .)  *XL“*i  /'  * 1 
GOT  06’  * 

MM  XL-  1*1  • 

XL0  w 1= 

GOTO*’  F 

6«4  <LX1  = ( '1  ,*T»)  ♦ S0»T  (T«  1)  t (T“l  + '911 

XL®  11=*!  *1*  (T7-T111  /(  2.  *7**1*  (TM*”!] 

6?*:  TF(  I*0M*  1,0*.  1 1 Of)  *’  •-?«, 

CY*  (XI  ■cot >4*  ’>  /?.-xtl  •- ? (Ml  - Y«*  r0(’l-l)  /’ , ) n* 

* Y=  (vi*'Kl*ll  /2.-Y10  »*3  (91-Yft»  TO  (1-11  /?. 1 /0Y 

nwi  = < V0*r3n*n  /? . -Y  1’  • S'*  1- Y »«  50 ( 9- l> /2 . ! /0Y 

•tyii*  (yi*;tk  i*  n / ?,  -th«£tii.  y3»ito(i-  n / ?.  t 7 iy 

00  TO  ft’’ 

626  'Y»  (Y1  *'■’('4*1)  /’.-YIO*®3  (Ml-  Y1»F3  (9-1)  /2.  t OY 

TYs  (Y9*T’(9M)  /’.-YU  **»  (11- Ya*TO(9-l)  /2.)  /9Y 
’Y^IM  y3  **»  Cl*  1)  /’.-Y  l'’*ril-Y3*EP(1-l)  /2.  1 79Y 
TTY  *n=  <Y  ?<►  £ T =»(•!♦  1 »/? . -*t9*ET91-  YS*£TO(‘l-l) /2. ) /TY 
62*  TFtniFT.i3.lt  00  TO  700’ 

CM2»Y’  'ftl(Y)  -Y  3*c‘ICJt 
-*?=Y’**n(9)-v  T»*l  (1) 

GO  TO  '=10 
750’  r**2  s’.*tocm> 

T9’  *2. **0(1) 

7505  rOMTTMi'T 

11* t = Y«»XL*  (2.  •<U'4  1‘  c“l/0Y-(  XL  41*  rYMl+r^iML’Il'TY-VIl)  1 
Hi  *XL*  XL4t*m  ‘Y7/OY*  »XL*(XL41*"741*E4t»X.941*Tf-Y41) 

1 Y1 0 *’  .*X“r  *OX  ’*  e,*l*  $ E®  + "°*  ( »Y1*  r«l  -T'12)  *X  1 

41(9,3'ftXLM?.*XL''l,:'<l»v6/0Y*(XL'*l»iY41*E91*XL’41*TY-V41>*Y9> 
11(9,1  )=-XL»E91*XLP41  *rv»Y? 

"1(9,  ?’s  OX  ’*X9  F.  ?,  •XL  •*>**•  YL  3U1  *rV*Y  10 
•1(9,3  It  XL*  Z U *XL’*U*  ”•  Y9 
01(9,1  1*01(9,31=0. 

’KM.’  ts-iy?*FY 

4 2(9,1  lx-’.  *XL*»AL»XL  Mt»5'il»TY»YS 

4 2(9,2  )»-(4.*XL»»AL»Xl'*l*£*U*FY»Y19*SEP*X*(Yl*T'*l-r921  > 

12(1,7)’’.  *XL»xiL,XL“l‘"'*l*rY»YQ 

•2(1.1  ) = XL*Y  ••  (’.•XL**  l*i*11/  <P°»9Y)  - (XLMt»FTY11*2,»XL*1l*STll*TV 
1 -»’*Y“1 1 /»»! 

•2(9,2'  = -(«.*YL*XL4t*"***l*Y7/<P9*0Y) *(XL11»TTYlt*2.*XLF11»iTll»TY 
1 J*  V “1 1 • X i,  * Y 10  *’.0/  OJ»0£P»X»»1»T11) 

•2(1,1  ) = XLM2,  •XL'*l*F'il*Y5/0Y*  (XL11*ETY'*1»2.  • X L"4 1**T  *«l»T  Y-PR* 

1 Xll)»Y11/®7 
02(9,2 )*-0x’*Tv 
02(9, 11=02(9. 11=0. 

41(1,1  )»13  (1,  3 )»B. 

4Y(9,2t*1X2*X*Yl 

53(1,1 1 = 11(1,2  1 = 33(1, 1’  *0, 

03(1,1  I*-XL*X9 
C3(1,2  )*-’•  *yL  *Y  10 
03(1,3  1*  XL*  Y9 

31(1)  «"V’»FY*  ( F11»XLB1t  •TY-V11)  «’ 1 1 * *2 • ( X9E *3X  2 ,X • X 1)  *SE° 
32(1>*',X’*(XL°',l*iTHl  *’Y/PP-  Vt*l  1 • TY  + 0X  2*  XAL  *Xl  **  1*  ’ 41*r  Y *»2-X*  Y 1 
l *T91*  r11*Orp 
33(  91=  **=92 
*014  C 09711  "F 
C 

C THE  FT  *1. 1 POINTS  ;VAL  l*  TION 

C 

c 

0 TME  0U’F9  -OSE  50U134  ?v  ’OIOITION 


005700 
005’10 
005720 
005*30 
005740 
005750 
005760 
005770 
005790 
005790 
005900 
005310 
005320 
005630 
005340 
005950 
005360 
005370 
005360 
005390 
005901 
005910 
005920 
005930 
005940 
005950 
005960 
005970 
005960 
905990 
006000 
006010 
006020 
006030 
006040 
* 006050 
006060 
006370 
006080 
006090 
006100 
006110 
006120 
005130 
006140 
006150 
106160 
006170 
006130 
006190 
106200 
006210 
006220 
006230 
006240 
006250 
006260 
006270 
006230 
006290 
006300 
006310 
006320 
006330 
006340 
006350 
006360 
006370 
006160 
006390 
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.*> 

106411 

V > 

**0  * *>  j ? •* ; f * 

106411 

Ml* 

ni  :rt  ** .:» * '- ■* « t * r » i>  * « it zn*, n * n ft-?-,*,  t > *■»»  t i-  n* ,n  *?3 < 

116421 

i tr *>'»*  .r<  *-i  <t  i)  ***  <:*i5c,  :>*:?  (fioe,  n *o. 

106430 

' 1 f '•'.Mi  1 . 5 

016<41 

006451 

■'K!")'’*)*!.  1 

116460 

?»(  t?5  *•?»■!,  1 

105470 

TF(5-*9».  ) 5 n T->  9"i', 

106491 

t' ?G€-i ) ) 

006490 

006500 

T*C-'I**.T?.t)  *“?**.  **0  (1*15*) 

006510 

13(1*1*-, '>  *1X?*<*»1 

006521 

m i*i**  ,i»  *’.*?»<••?*»./(  l,  ♦**<> 

006531 

<13(  15?  rs  ,?)*-?,»  VX<»  ( 1.  »VX<)  • *(. 

006541 

"■*f  * ?•  •*  •(  ? . *XX  <f  1 .)  / ( 1 .4  xx  <> 

006551 

13( 

006560 

-o  Tn  a*"4* 

916571 

VR1  *73  <T  > 

006591 

Tci it3 *,T i « ' t . i ) v*n»v9ir"or,-> 

005591 

“3(T*1  f'.'iu.  5 

006601 

•y»( 

056610 

on* 

IX-TV't" 

005620 

r 

006631 

c 

tut  OMT’  in'-.-  90U5I3»9»  50NOITION 

006640 

’6T'J9'4 

006650 

:H> 

006660 

SW»"II t*I*  Of»HrMS(T»niFi,:p,ir;,X5,-<5TS5T,i:N,H,S,r,  9LT) 

016671 

1<X»x<W  r,,  39,  or v , KMT-IC,  01*5 s,  35,  IX,  cm,  M'S  *11,  71515,  *i 

, 106591 

i pf  ,"*  ,u*,7i$fiF,S'j,r »*,5s,3»w,  5i,3  9909,'Ctr4, i'i5:,r;'ni, ivraer, 

116490 

? o?T,y  v<,3T9y,  »u  11,71  "3  6 2 , *1  IT*",  5*°3,  JUS  ( J)  , T^l  ( 3)  , £ 3 ( 2 5 3 > , 

006711 

3 *N  (71  9)  ,*?  (?0  0)  ,-r rj  ( ??1)  ,rr*|{?11)  , **»  (?10>  ,*1  ( ’"O)  ,*1  ( 219)  , J211, 

016711 

4 F»(?1  «>  ,’9(210)  ,13(2  0")  ,X»lr((211)  ,VN<?"9>  ,73(211)  , 7»(  31  0 ) , *»  ( 21 ) > 

, 106721 

1 

5 01  (’1  "1  ,"?(?1  9)  ,53(?  35) 

005730 

It  i 

ITi'lS*"'!  7 (?1">  ,Z  (7  , 16  1 

006740 

25 

**>4.1’  ( 1X5,  .5  X,  I'MO’TIL'  *1*  5T4TJX  $ **14.8) 

156750 

40 

'OJXS*  (1  X9«|*  15**  .4  1 

006760 

4i 

£Ti*  15P*.-  > 

106775 

4? 

*0R4ST  (••*  fir  15*1  ,4  1 

006*91 

43 

*0»11T  (94  T 1 * 15*9  ,4  ) 

006731 

66 

*0*4»T  (99  71*  15*9 .»  1 

506190 

46 

*0»9»r  f* M *3*  13*1  ,3) 

006910 

517 

*0’«MT(1'<  Y/HT«  15*»,4  1 

004920 

5 5» 

*0»91T(04  ?0/3 0**15* • ) 

006930 

510 

*034ST  (*4  «4*m*  15*1, ♦) 

906441 

511 

*0»«X4T  (•  4 9T/ooo*i5*«  1 

J0605  0 

51? 

•fr»54r(«'X  9T/?r*  n**,;  > 

506861 

51? 

*0341*  (5U  4/w*  * 15*9. x 1 

006070 
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Appendix  C 

Four  Key  Subsystems  Within  I tract 


Nondimensionalizing  the  Variables  and  Initializing  the  Grid 

Prior  to  entering  the  computational  loop  the  working  variables  were 
nondimensionalized  or  normalized.  These  variables  were  listed  below 
along  with  a definition  of  each.  The  format  selected  was  to  present  the 
coded  variable  on  the  left  side  of  the  equal  sign  and  the  real  or 
physical  definition  on  the  right  side  of  the  equal  sign.  No  explanation 
was  Included  as  to  choice  of  normalizing  factors. 
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With  Eq  (54)  defined  for  all  cases,  some  others  depended  on  the  value  of 
w.  If  a)  were  not  equal  to  zero,  then 
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where  the  reference  temperature  was  taken  as  TW(Y-1 However,  for 
the  case  where  w was  equal  to  zero,  the  quantities  of  Eq  (55)  plus  one 
were  defined  as  follows: 
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These  quantities  were  frequently  used  in  the  grid  computation  and  pro- 
vided a summary  of  the  nondimensionalizing  techniques  used  throughout 
the  code.  Before  beginning  this  computation  within  the  grid,  however, 
there  had  to  be  an  Initialization  of  the  profile. 

Initialization  began  by  defining  Y in  the  code  as  the  distance 


AnK+i) j-1 


An-j  which 


measured  along  the  n axis.  Any  An*  was  defined  as  , . 

3 [ AnK 

yielded  a fine  mesh  of  nodal  points  near  the  surface  and  an  adequate 
spacing  toward  the  edge.  Y values  were  assigned  by  successively  adding 
all  An  values  from  the  surface,  to  the  point  In  question.  Then,  three 
hypothetical  successive  columns  of  nodes  were  created  by  the  following 
statements: 

D1  s D2  = D3  = 0.,  from  the  surface  to  the  edge  of  the  boundary  layer. 
Incorporating  the  notation  of  fig  1, 

Vj  j * Vl-l,  s V1-2,  " “Yj*  *or  J from  surface  to  the  edge  of 
j J 

the  boundary  layer. 

In  a similar  manner,  three  successive  stations  of  F,  6,  <e,  and  i were 
assigned  values  of  1.0.  Finally,  all  coefficients  of  the  system  of 
finite  difference  equations  were  set  equal  to  0. 

This  Initialization  provided  the  primer  to  begin  the  backward 
differencing  along  the  £ direction  and  the  central  dlffefenclng  along 
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the  n direction.  The  finite  differencing  system  was  unconditionally 
stable  for  Increments  of  An  and  A£,  and  the  iterative  stepping  procedure 
along  £ damped  the  error  due  to  the  grid  Initialization  within  a few 
steps  (Ref  6). 


The  Finite  Difference  System 

Coefficients  of  the  finite  difference  equations  were  computed  for 
the  matrix  equations  which  would  be  solved  In  a succeeding  step.  These 
equations  were  derived  starting  with  the  concept  of  a grid  as  in  fig  1, 
and  the  stipulation  that  a function  could  be  described  at  a point  by  a 
Taylor  series  expansion  about  another  point.  For  Itract  the  approxima- 
tion was  made  that  for  any  functional  value,  F, 
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the  second  and  first  partial  derivatives  of  Eq  (57)  were  expressed  by 
central  differencing  as  follows: 

32F(1,j)  _ Y6F(T,j+U  . 2Y7F(i,j).  + Y8F(1J-Jl 


3F(1,J)  . Y9F(1,J+11  . Y 
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3V  30 

The  same  format  of  expression  was  used  for  ^ and  — -y.  For  a 


streamwlse  series  of  nodal  points  along  e the  backward  differencing 
system  was  written  from 
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F(1-2.J)  - F(1  ,J)  - (A?,. §+ 21  ^ 

Only  expressions  for  the  first  derivative  with  respect  to  5 were  required 
and  this  equation  was  as  follows: 
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36 

Again  the  same  format  of  expression  was  used  for  and  all  derivative 
forms  of  Eqs  (20),  (21),  and  (22)  had  been  derived.  Then,  due  to  their 
recurring  use,  the  following  definitions  were  made  for  computational  con 
venience  and  efficiency: 
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(62) 


Through  Taylor  series  expansions  about  F(1,j)  and  T(i,j)  and  neglecting 
terms  with  second  order  partial  derivatives  and  higher,  then  FM1  and  TM1 
were  actually  expressions  for  F(1,j)  and  T(1,j),  respectively. 

Returning  to  Eqs  (20),  (21),  and  (22)  and  the  construction  of 
linearized  finite  difference  equations,  there  were  three  types  of  non- 
linear terms  with  which  to  be  dealt.  Using  F and  G to  represent  any  two 
general  function  symbols  the  nonlinear  terms  were  of  the  types: 

(F)(||],  and  (F) (G) . where  F could  have  been  equal  to  G. 

Returning  to  the  notation  of  the  problem  variables  It  was  shown  that 
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where  FY  was  equal  to  ^ a known,  and  — was  unknown,  and 

that 


r - 2F(1,j)  f(1-l,j)  - F(1-l,j) 


(65) 


where  only  F(1,j)  was  unknown.  All  terms  had  been  represented  in  finite 
difference  form,  and  the  final  step  incorporated  these  linearized  models 
Into  Eqs  (20),  (21),  and  (22)  to  derive  the  overall  system  of  finite 
difference  equations  (Ref  8:67-71). 

From  this  system,  the  coefficients  of  F(i,j-1),  F(i,j),  F(1,j+1), 
and  T and  V at  these  stations  were  collected,  computed,  and  passed  to 
the  matrix  Inversion  routine  resulting  In  solutions  for  F,  V,  and  6 from 
the  surface  to  the  edge  of  the  boundary  layer  at  the  current  station,  s^. 


Subroutine  Reystr 

This  routine  was  called  from  the  main  program  at  each  station,  s^ , 
at  and  beyond  the  point  of  transition  to  turbulence.  The  purpose  of  this 
subroutine  was  to  calculate  an  eddy  viscosity  for  the  inner  and  outer 
regions  of  the  two-layer  turbulent  boundary  layer  model. 

Computation  within  Reystr  began  with  Taylor  series  expansions  of  F 
to  the  third  order  partial  term  about  the  first  station  at  the  wall. 

With  values  for  Fj^,  F^.  ^=3  an(1  ^j=4»  8 four-point  finite  difference 
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, and  the  coefficients  of  the  F terms  at 


3F 

expression  was  formed  for  ^ 
each  node,  one  through  four,  were  represented  by  Yll,  Y12,  Y13,  and  Y14 
In  the  code.  Next,  a nondlmenslonal  molecular  viscosity-density  term 
was  calculated  for  the  wall  with  a shear  stress  term  that  followed: 
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PI2  - 3F 

PI2 ' TwTTst;*, 


An  iterative  loop  was  begun  to  generate  the  nondimensionalized  inner 


eddy  viscosity  model. 


'inner 


, of  Cebeci-Smith-Mosinskis  for  each  node  in 


the  n direction  for  the  current  s^.  In  the  actual  code  and  following  the 
calculation  of  a number  of  Interim  quantities  that  did  not  necessarily 
represent  any  real  boundary  layer  characteristic,  three  important  computa- 
tions were  made.  First,  6/L  was  calculated.  Next,  an  intermediate 
quantity,  DD,  to  be  used  later  in  the  outer  eddy  model,  was  calculated. 
Finally,  PI1,  another  Intermediate  quantity  used  in  the  inner  model,  was 


computed: 


6/L  = XNN. 


U-  .995llAXKNJ.(J.1)l 


j-U-l) 


edge  of  the 
boundary  layer 
DD  = l 

j=2 


e e 


. “bilfLhil 


[(y-Dh2]1-5 


T_(y-DH;+198.6. 


P7F- 

i 5 


where  the  s was  Included  for  the  case  of  conical  flow  only.  Again,  a 

jr 

term  was  generated,  but  using  only  a three-point  central  differencing 


92 
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scheme  on  this  occasion.  The  final  step  of  the  loop  was  the  actual  com- 
^ i nner 

putation  of  — at  the  current  node  j: 


6inner| 

Pj 


. 1 6 ( P II ) (1  -exp(-[(PIl)(PI2)]1/2/(26 


(PP) 


(py) 


j-)))1 


(69) 


™ F(1>j*1}  - Y1°  F<i»j)  - T 


An  4 


f(pp)i 


where  Y8,  Y9,  and  Y10  were  coefficients  obtained  through  Taylor  series 
expansions  of  F(i,j-1)  and  F(i,j+1)  about  a point  F(i,j).  As  the  calcula- 


e 


tlon  of 
e1nner 


inner 


progressed  from  the  wall  out  into  the  field  of  flow. 


was  greater. 


« • 

j+1  retained  its  own  computed  value  or  that  of  1nnerj , whichever 


The  outer  law, 


outer 


, was  computed  through  an  iterative  loop 


similar  to  that  of  the  inner  model.  It  culminated  with  the  expression 


e°uter  _ 

2XRe_ 

00 

1/2 

DD 

^ . U 1 00 

Ta>+198.6 

ye 

yref 

J 

(py)j 

Twi; 

Ii 

lTeJ 

2 

s 

UY-l)Ml+198.6 

. - 

(70) 

where  the  s was  Included  only  for  the  case  of  conical  flow.  In  order  that 
a compatible  combination  of  computed  viscosities  were  retained,  the  values 
of  eddy  viscosity  from  the  outer  law  replaced  th  2 of  the  Inner  law  from 
the  point  of  Intersection  of  the  graphs  to  the  edge  of  the  boundary  layer. 
Graphically,  this  was  depicted  In  Fig  18. 


t 
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Fig.  18.  Matching  the  Inner  and  Outer  Eddy 
Viscosity  Models  (From  Ref  8:21) 


Having  calculated  the  initial  eddy  values  for  the  inner  and  outer 
viscous  regions  of  the  boundary  layer,  it  was  appropriate  to  subject  this 
model  to  two  more  factors.  Both  were  factors  of  degradation  and  were 
included  to  better  describe  the  character  of  turbulent  activity  within  the 
boundary  layer. 

Objections  have  been  raised  to  the  use  of  an  eddy  viscosity  term,  e, 
in  place  of,  or  in  addition  to  the  molecular  viscosity,  y,  of  a fluid,  u 
is  a real  property  of  a fluid,  e is  only  an  effective  description  when  a 
fluid  Is  in  motion,  and  It  Is  clearly  not  a property  of  the  fluid.  But, 
with  reservation,  it  has  been  used  to  express  the  behavior  of  turbulent 
stresses  in  terms  of  mean  velocity  gradients  of  a flowing  fluid.  It  has 
been  possible  to  obtain  a satisfactory  description  of  mean  properties 
within  turbulent  flows  by  assuming  this  flow  to  behave  as  a Newtonian 
fluid.  Incorporating  an  eddy  viscosity  model  along  with  y,  and  including 
two  factors  of  Intermlttency  when  appropriate  (Ref  20:25-26).  A laminar 
and  Irrotatlonal  flow  became  turbulent  as  It  passed  through  a region  of 
transition  In  which  only  a fraction  of  the  time  was  spent  In  a turbulent 
state.  During  that  time  In  laminar  motion,  the  Reynolds  stress,  hence  e, 
would  have  been  zero.  Then,  to  adequately  describe  the  effects  of  e at 
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any  point  by  the  relative  fraction  of  time  that  that  point  would  be 
engulfed  in  turbulent  flow  (Ref  21:117).  Therefore,  the  first  multiplica- 
tive factor,  called  an  intermittency  factor,  was  applied  to  e to  more 
accurately  describe  the  e within  the  transition  region.  The  intermittency 
or  probability  factor  of  Dhawan  and  Narasimha  was  used  for  this  program. 
The  factor  was  computed  as  follows  (Ref  8:28-29): 


r(s)  = 

1 - exp 

* 

-.412 

^ 1 

» 

Then,  the  computed  ^ 


scurrent  ” transition  point 

r— p.  C 

v"  'transition  point  J 
original  was  replaced  ** 


"2l 


(71) 


modified 


= (T(s))  ^ , . , 
''  P|Onginal 


(72) 


The  second  factor  was  then  considered.  It  was  observed  by  Klebanoff 
that  in  a turbulent  boundary  layer  with  a free  boundary,  as  the  free 
stream  was  approached  the  turbulence  became  intermittent.  This  inter- 
mittent nature  was  observed  first  at  y/6  greater  than  .4  with  less 
turbulent  intensity  as  y/6  grew  larger.  It  was  thought  that  a good 
prediction  of  turbulent  intensity  probably  depended  on  a correct  weighting 
of  the  probability  density  for  the  turbulence  of  the  free  stream  with  that 
within  the  boundary.  It  was  found  that  a good  description  of  Y'  was  a 
Gaussian  Integral  curve  given  by 

Y'  - 7 (1  -erf(C'))  (73) 

where 

«'  ■ (•*  I)"’  [J-  ■ 5[f-  .78]  <7«) 


M. 


95 


These  expressions  Indicated  that  the  edge  of  the  boundary  layer  had  a 
random  character  with  a mean  position  at  y/6  equal  to  .78.  The  edge 
vacillated  from  y/6  equal  to  .4  to  y/6  equal  to  1.2.  Finally,  If  It  were 
assumed  that  the  free  stream  contributed  little  to  the  measured  turbulent 
quantities  of  the  boundary  layer,  an  allowance  could  be  made  for  the 
effect  of  Intermlttency  by  dividing  by  the  factor  if'  (Ref  22:15-18). 

Cebecl  used  the  approximate  expression  for  Eq  (73)  to  give  a multi- 


plicative version: 


1 +5.5 


(Ref  7:1679) 


which  led  to  the  coding  for  this  second  factor.  If  y' were  not  included, 
then  a newly  defined  viscosity  was 


e = 1 + - r(s) 

u 1 


Including  y'»  Shang  formed  the  following  model: 


e * 1 + 


1+5.5 


+ 1 i 


2*75 


For  purposes  of  this  study  Eq  (76)  became  eddy  model  zero,  and  Eq  (77) 
became  eddy  model  one.  Then,  whether  or  not  y*  was  Included,  the  quantity 


e was  defined  by 


• ■ 1 ♦ pjr-  (e  - 1) 


In  a final  note,  the  decision  of  whether  to  use  eddy  model  zero  or  eddy 
model  one  depended  on  the  original  assumption  that  either  the  free  stream 
turbulence  had  an  effect  on  the  e of  the  boundary  layer,  or  It  did  not. 
This  factor,  y'»  was  to  have  a definite  effect  on  the  analytical  results. 


•i 


r>->  «.*;> 


and  this  entire  subroutine  was  Included  with  the  program  listing  of 
Appendix  B. 


Subroutine  Cfstno 

Like  Reystr  this  routine  was  called  from  the  main  program.  But 
unlike  Reystr,  Cfstno  performed  Its  computation  throughout  the  laminar, 
transition,  and  turbulent  regions  of  flow.  The  purpose  of  this  routine 
was  to  calculate  a Stanton  number,  a measure  of  heat  transfer;  the  local 
coefficient  of  friction.  Indicative  of  shear  stress  at  the  surface;  and 
Reynolds  numbers  based  on  displacement  thickness  and  momentum  thickness. 

(pm)w 

Computation  began  with  yp-^y-,  coded  XLM1  in  the  program.  The  formula 
by  which  XLM1  was  computed  depended  on  the  value  of  the  exponent  in  the 
viscosity  law  of  Sutherland,  the  value  of  this  exponent  being  specified 
by  the  programmer.  If  the  exponent  were  zero,  then 


XLM1  = ==• 


Tjl/2  fT  +198.6 


If  this  exponent  were  one,  then  XLM1  was  one.  Otherwise, 


T w-l 

XLM1  « J-  (80) 

. e 

Next  to  be  calculated  were  transformed  quantities  similar  to  4 or  heat 

flux  and  T or  shear  stress.  First,  the  same  four-point  finite  differ- 

9FI 

ence  scheme  used  In  Reystr  for  was  repeated  at  this  point  to  calcu- 
late |^|w  and  f^|w*  Then  the  transformed  t,  coded  TAU,  was  computed: 

(pvL  P«  ue  fuj2  aFi  1 

q (81> 


For  the  case  of  the  axl symmetric  flow,  both  TAU  and  QS  were  divided  by 
the  nondimensional  station,  s..  With  this,  preliminary  calculations 
were  completed. 


A Stanton  number  and  coefficient  of  friction  followed  next  In  the 
computation.  If  Tw  equaled  T0,  there  was  no  heat  transfer  and  St,  coded 
STNO,  was  zero.  Otherwise, 


The  model  from  which  this  expression  came  was 


e P„u0(H  -ti  ) 
e e e w 

For  the  calculation  of  cf|local  stat1on.  coded  CFNO, 


^refl1/2M.Ll-1/2  <^w 


PJVef  u» 


With  St  and  cflocal  computed,  only  the  transformed  expressions  for  Reg* 
and  Ree  remained.  Coded  as  REYDT  and  REYMT,  these  quantities  were  com- 
puted from  the  following  statments: 
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This  completed  calculations  within  this  routine,  and  further,  com- 
pleted the  formal  description  of  four  important  subsystems  within  Itract. 
Again,  this  subroutine  was  Included  with  the  program  listing  of  Appendix 
B.  In  this  appendix  consideration  was  given  to  the  important  concepts  of 
the  nondimensionalizatlon  of  working  quantities.  Initialization  of  the 

grid,  and  the  generation  of  finite  difference  coefficients.  Also 

• 

Included  was  a brief  description  of  the  two  subroutines  used  in  the  compu 
tation  of  eddy  viscosity,  heat  transfer,  and  skin  friction.  The  theory 
presented  In  this  appendix  should  provide  a better  understanding  of  the 
code  In  general,  and  the  modification  for  mass  transfer  specifically. 
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Appendix  D 


Fortran  Computer  Code  Key 


Coded  Symbol  Represented  Quantity  (Values  included  for  those  quantities 
remaining  constant  throughout  this  project) 


Inputs 

(In  order  read  by  computer) 

G 

Y = 1.4 

PR 

Pr  = .73 

XMINF 

TA 

T„ 

DS 

Stepping  Increment  In  s,  along  the  streamwise  direction, 
DS  = .0004 

SI 

Initial  station,  s, , began  computation  within  the  grid, 
SI  = .0006  1 

OMEGA 

Exponent  in  the  viscosity  law  of  Sutherland,  OMEGA  = 0 

ERROR 

A convergence  criterion,  the  acceptable  difference  between 
the  quantity  w calculated  In  two  successive  calls  of 

the  matrix  Inversion  routine  at  the  same  station  s^ 

XXK 

^i+l 

— , a constant  ratio  from  surface  to  the  edge  of  the 
flnj  boundary  layer 

BO 

Tw 

BTRX 

Station  s,  at  which  transition  from  laminar  to  turbulent 
flow  began 

• 

PRT 

Prt  * .9  (exceptions  noted) 

XINTER 

A flagged  quantity;  XINTER  = 0.,  eddy  model  zero  was 

used;  XINTER  ■ 1.,  eddy  model  one  was  used  In  the  compu- 
tation of  e 

DYW 

An^,  the  first  Increment  In  n 

I EDGE 

Total  number  of  nodal  points  or  divisions  In  the  n direc- 
tion within  the  grid 

3 

INTACT 

Not  used  In  this  study 

100 


IDIFF 

I 


IEND1 


MSP 


J2DA 


IPRES 


I CHS 


IPRN 

I 

XLGTHMD 

RINFA 

IBLW 

STRT,  DONE, 
RVRAT 


NUMDAT, 

XPOS, 

RHOVRAT 


A flagged  quantity;  IDIFF  = 0,  a three-point  differencing 
scheme  was  to  be  used;  IDIFF  = 1,  a two-point  differen- 
cing scheme  was  to  be  used;  IDIFF  was  set  equal  to  0 for 
this  project. 

Total  number  of  nodal  points  or  divisions  in  the  £ direc- 
tion within  the  grid 

A flagged  quantity;  MSP  = 1,  program  printed  abbreviated 
data  from  each  station  computed;  MSP  = 5,  program  printed 
every  fifth  station;  MSP  was  set  equal  to  1 for  this  pro- 
ject. 

A flagged  quantity;  J2DA  = 0,  designated  a flat  plate 
calculation,  J2DA  = 1,  designated  an  axisymmetric  cone 
calculation 

A flagged  quantity;  IPRES  = 0,  indicated  that  dp/dx  was 
zero;  IPRES  = 1,  indicated  that  dp/dx  was  not  zero;  IPRES 
was  set  equal  to  zero  for  this  project. 

An  array  of  integers  which  designated  stations  where  a 
double  step  was  to  be  taken  between  computations  of  a 
column  of  nodal  points 

An  array  of  integers  which  designated  stations  where  a 
full  profile  of  boundary  layer  data  was  to  be  printed 

Length  of  the  model , L 


P. 

A flagged  quantity;  IBLW  = 0,  no  mass  transfer  consid- 
ered; IBLW  = -1,  mass  transferred  at  a constant  rate; 
IBLW  = 1,  mass  transfer  varied  along  the  length  of  the 
model 


If  IBLW  * -1 , mass  transfer  began  at  some  number  of  feet 
from  the  leading  edge  or  tip  and  continued  to  some  other 
position  downstream,  transferring  at  a constant  rate, 
(pv)  (pv) 

(■p-jyw-  for  the  plate  or  for  the  cone 
00  0 


If  IBLW  = 1,  this  stipulated  a varying  transfer  rate 


pv. 


pu„ 


posl 


or 


beginning  at  X , at  a strength  of 

(pv)w  ^ _ . 

and  continuing  to  X ..  - at  a correspond- 

(Pue)p0si  P°s  f1nal 

Ingly  specified  strength.  Varying  transfer  rates  were 


designated  In  between. 


101 


Outputs  and  Miscellaneous  Working  Quantities  Applicable  in  This  Study 

(alphabetical ) 


BCVW 

BLDT 

BLMT 

BLT 

CFNO 

EO 

ETA 

FI 

H/HE 
MACH 
N,  XNN 
RE 
REY 

REYEXT 

REYDT 

REYMT 

RO/ROE 

RVRAT(VRVRAT) 

STNO 

T1 

TRR 

TW/TE 

UE 


V(1 ,1 ) , defined  In  Eq  (29) 

6*/L 

0/L 

5/L 


Cflocal 

Eddy  Viscosity,  either  from  eddy  model  zero  or  eddy 
model  one 

n 

u/ue  (in  the  output  listing  only) 

H/He  (In  the  output  listing  only) 

Mach  number  (in  the  output  listing  only) 
y/L  (listed  as  N in  the  output) 

P/Pe 

Re 

00 


Re6* 

Re0 

P/Pe  (in  the  output  listing  only) 


(pv) 


w 


(pv). 


current  station 


or 


T*rr 


Tpsi; 

Ste 

T/Te  (In  the  output  listing  only) 
(V198.6)/(T„(Y-l)Mf+198.6) 


current  station 


VT.  (In  the  output  listing  only) 
u/ue 


VI 


Y/BLT 


v/ue  (In  the  output  listing  only) 
Defined  in  Eq  (27) 

3 

station  s^ 


y/6  (in  the  output  listing  only) 


Appendix  E 

A Cubic  Spline  Approximation  for  the 
Description  of  General ly  Varying  Mass 
Transfer  hate 

In  modeling  or  mathematically  describing  a varying  mass  transfer 
rate  It  was  assumed  that  through  some  means,  there  would  be  knowledge 
of  the  strength  of  mass  transfer  at  a finite  number  of  stations  along 
the  model.  So,  there  was  information  of  the  form  (x^,f^)  for  1 values 
from  1 to  n.  The  objective  was  to  construct  a function,  f(x),  such  that 
fj  was  equal  to  f(x^)  and  that  f(x)  was  twice  differentiable  over  [x^ ,xn] 
This  f(x)  would  provide  the  value  of  mass  transfer  for  any  station,  s^, 
along  the  surface  of  the  model.  Figure  19  depicted  the  curve  to  be 
specified. 


Fig.  19.  Building  a Cubic  Polynomial 
Between  Any  x^  and  x^+^ 

The  function,  f(x),  was  specified  as  a different  cubic  polynomial 
In  each  Interval,  x^  to  x^+j.  It  was  required  that  the  function  be  con- 
tinuous, together  with  Its  first  two  derivatives,  at  each  Junction 
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between  two  polynomials.  Thus,  for  each  [x^.x^],  an  f(x)  was  constructed 

equal  to  l C1^  . The  function  was  formed  recursively.  Supposing  that 
1=0  *- 

f(x)  had  been  generated  to  an  x equal  to  x- , It  was  necessary  to  choose 
C^_q  j such  that  f(x),  f'(x),  and  f"(x)  were  continuous  at  x^ , and  it 
was  left  to  find  f(x)  over  the  interval  [x1 »x1+1 ] . This  led  to  four  linear 
algebraic  equations  with  four  unknowns,  tQ  3.  These  equations  were 

as  follows:  , 


f(x4) 


C1 

xi 


f(*i>  ■ c<  ’ 

f"(x.)  « l cl£(£-l)xf- 

1 £=2  L 1 


f = f C1 
1+1  1 1+1 

This  system  was  solved  up  to  x^+-|  at  which  point  the  process  was  repeated 
from  xi+1  to  x1+2  (Ref  11). 

Returning  to  the  Initiation  of  this  recursive  procedure,  values  were 


known  for  (x4,fj)  for  1 equal  from  1 to  n.  Then,  fi  was  approximated  by 
f2"fl 

— — and  fV  was  approximated  by  the  expression  ■ ■ - . The  Initial  condl- 

x2  X1  x2”xl 

tlons  were  then  f-j , fj,  fj,  and  fg.  The  four  equations  Initially  to  be 

solved  were,  then,  given  by 


C0  + C1  X1  + C2  xf  + C3  x? 

c]  + 2Cg  x1  + 3cJ  x? 


24  + 6C3x1  “fl 
ej  ♦ c]  x2  ♦ cj  4 ♦ c’  4 . f2 
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(89) 


Solving  for  c]  yielded  a cubic  polynomial  expression 

fj(x)  = cj  + c]  x +•  C J x?  + C3  x3 

which  was  descriptive  of  an  appropriate  curve  connecting  points  one  and 
two.  Then,  having  specified  the  polynomial  for  this  first  interval,  the 
successive  polynomials  and  their  intervals  were  recursively  computed  to 
xn  as  previously  discussed,  though  now  a polynomial  expression  existed 
for  finding  f'(x)  and  f"(x). 

Finally,  then,  for  any  position,  s',  along  the  surface  of  the  model, 
the  Interval  s^  to  s1+l  1n  which  the  position  was  contained  could  be  found. 
Knowing  the  interval  was  to  also  knew  the  corresponding  cubic  polynomial 
that  described  that  increment,  and  hence,  the  value  of  mass  transfer  rate, 


A Heat  Transfer  Comparison  with  Moffat 


and  Kays,  Mass  Transfer  Factor  of  700115 


Table  VIII 


A Heat  Transfer  Comparison  with  Moffat 
and  Kays  at  the  Suction  Asymptotic  Limit 


Itract  Prediction, 
Stoo(lO)3 


stoker 

Experi- 

mental 

Eddy 

Model 

One 

Eddy 

Model 

Zero 

9.33 

4.90 

5.00 

8.07 

8.34 

8.38 

7.75 

8.17 

8.20 

7.82 

8.09 

8.12 

7.64 

8.04 

8.06 

7.99 

8.00 

8.03 

7.71 

7.98 

8.00 

7.85 

7.96 

7.98 

7.82 

7.95 

7.96 

7.95 

7.94 

7.96 

7.94 

7.93 

7.95 

7.91 

7.93 

7.94 

8.24 

7.92 

7.93 

8.17 

7.92 

7.93 

7.82 

7.92 

7.92 

7.97 

7.91 

7.92 

7.88 

7.91 

7.92 

8.35 

7.91 

7.92 

7.76 

7.91 

7.92 

7.97 

7.91 

7.91 

7.75 

7.91 

7.91 

7.75 

7.91 

7.91 

8.08 

7.91 

7.91 

7.85 

.... 

.... 

Percent 

Error 

with 

Eddy 

Model 

One 


(Ref  13 
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Table  IX 


A Heat  Transfer  Comparison  with  Moffat 
and  Kays,  Mass  Transfer  Factor  of  .0019 


Appendix  G 

Cone  Heat  Transfer  Data 


Table  X 

A Heat  Transfer  Comparison  with  Martellucci, 
Laganelli , and  Hahn,  Data  Group  132 


Station,  s 

St  no)4 

st  PeUeCiq>4 

epmu„  1 
Itract 

Percent 

Error 

with 

Eddy 

Model 

One 

rheoretical 

Fully 

Laminar 

Experi- 

mental 

Theoretical 

Fully 

Turbulent 

.191 

4.0 

4.77 

7.7 

4.63 

2.9 

.262 

3.5 

3.50 

7.3 

3.96 

11.6 

.315 

3.0 

2. 4,3.1 

7.0 

3.78 

18.0 

.399 

2.8 

3.93 

6.8 

5.33 

26.3 

.470 

2.6 

5.67 

6.5 

6.86 

17.3 

.542 

2.4 

5.68 

6.3 

7.55 

24.8 

.589 

2.3 

6.91,6.61 

6.0 

7.67 

9.9 

.607 

2.2 

6.18 

6.0 

7.70 

19.7 

.732 

2.0 

6.79,6.84 

6.23,7.30 

5.9 

7.47 

2.3 

.750 

2.0 

7.19 

5.9 

7.41 

3.0 

.816 

1.9 

7.67 

5.8 

7.23 

5.7 

.958 

1.8 

6.31 

5.6 

6.89 

8.4 

(Ref  14) 
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Table  XI 

A Heat  Transfer  Comparison  with  Martellucci, 
Laganelli,  and  Hahn,  Reference  Data  150 


Station,  s 

St  (10)4 
00 

Experimental 

Ste  ^OO)4 
Itract 

Percent 

Error 

with  Eddy 
Model 
One 

.191 

2.10 

2.06 

1.9 

.227 

2.70 

2.63 

2.6 

.263 

2.14 

3.76 

43. 

.317 

1.91 ,4.03,3.86 

5.00 

19.4 

.353 

6.48 

5.33 

17.7 

.400 

5.31 

5.42 

2.0 

.544 

4.90 

5.08 

3.5 

.592 

4.30 

4.97 

13.5 

.610 

3.94 

4.93 

20.1 

.645 

3.80 

4.87 

22.0 

.681 

3.53 

4.80 

26.4 

.717 

2.12 

4.74 

55.3 

.735 

4.27,3.98 

4.10,5.08 

4.73 

6.9 

.819 

4.19 

4.62 

9.3 

.890 

5.99 

4.54 

24.2 

.962 

4.69 

4.46 

4.9 

Table  XII 


A Heat  Transfer  Comparison  with  Kartell ucci, 
Laganelli,  and  Hahn,  Reference  Data  1 


Station,  s 

St  (10)4 
00 

Experimental 

Ste  j^pdo)4 

CO  00 

Itract 

Percent 

Error 

with 

Eddy 

Model 

One 

.173 

3.73 

3.97 

6.0 

.191 

3.85 

3.92 

1.8 

r*. 

CM 

C\J 

• 

4.58 

5.03 

8.9 

.263 

5.33 

7.03 

24.2 

.317 

8.84,7.90, 

7.18,7.53 

8.79 

.6 

.353 

9.68 

9.22 

4.8 

.400 

9.12 

9.30 

1.9 

.472 

8.15 

8.99 

9.3 

.544 

7.69 

8.64 

11.0 

.592 

7.19 

8.43 

14.7 

.610 

6.89 

8.36 

17.6 

.645 

7.49 

8.25 

9.2 

.681 

7.04 

8.13 

13.4 

.717 

6.99 

8.04 

13.1 

.735 

6.94,6.84, 

6.75,6.86 

7.99 

13.1 

.753 

6.84 

7.94 

13.9 

.819 

6.64 

7.80 

14.9 

.890 

6.79 

7.67 

11.5 

.962 

6.34 

7.54 

15.9 

(Ref  15) 
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Table  XIII 


i * 


A Heat  Transfer  Comparison  with  Martellucci, 
Laganelli,  and  Hahn,  Data  Group  60 


Station,  s 

stjio)4 

Experimental 

Ste  ^10)4 
I tract 

Percent 
Error 
with  Eddy 
Model 
One 

.226 

3.80 

5.57 

31.8 

.262 

3.09 

4.7? 

34.5 

.315 

2.45,2.49,2.54 

3.85 

34.0 

.399 

1.99 

3.77 

47.2 

.542 

1.78 

6.94 

74.3 

.589 

1.98 

7.32 

73.0 

.607 

2.17 

7.40 

70.7 

.648 

2.46 

7.48 

67.1 

.732 

1.87 

7.43 

74.8 

.750 

1.76 

7.40 

76.2 

.816 

2.00 

- 

7.28 

72.5 

.958 

2.31 

6.84 

66.2 

Table  XIV 

A Heat  Transfer  Comparison  with  Martel lucci, 
Laganelli,  and  Hahn,  Data  Group  203 


Station,  s 

St  (10)4 
00 

Experimental 

Ste  ’YOO)4 

^OO  oo 

Itract 

Percent 
Error 
with  Eddy 
Model 
One 

.263 

I.OO(IO)"2 

5.15O0)"1 

48.5 

.317 

3.42 

1.68 

50.9 

.353 

4.72 

2.16 

54.2 

.400 

5.13 

2.52 

50.9 

.472 

3.71 

2.60 

29.9 

.544 

4.02 

2.52 

37.3 

.592 

4.15 

2.44 

41.2 

.610 

3.74 

2.41 

35.6 

.645 

4.05 

2.35 

42.0 

.681 

3.65 

2.30 

37.0 

.735 

3.99,3.88 

3.98,6.41 

2.24 

42.3 

.753 

3.02 

2.21 

26.8 

.819 

3.96 

2.15 

45.7 

.890 

5.35 

2.09 

61.0 

.962 

4.33 

2.02 

53.3 

VITA 
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